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1 .  INTRODUCTION 


Overview 

The  response  of  multi-layered  laminates  can  be  obtained  from 
classical  laminate  theory  [1].  Once  the  mechanical  properties  of 
individual  unidirectional  laminae  are  known,  the  specific 
eigenvalue  equation  governing  the  frequency  response  of  the 
laminate  can  be  expressed  in  either  discrete  or  distributed 
parameter  form.  According  to  classical  laminate  theory,  the 
frequency  response  of  symmetric  composite  laminates  depends  upon 
the  stiffness  properties  (Du,  D12,  DM,  DM,  D,6,  and  D26)  .  Likewise, 
damping  can  be  accounted  for  by  appropriately  introducing  damping 
parameters  C^.  In  turn,  the  stiffness  parameters  D;j  and  damping 
parameters  are  known  functions  of  the  individual  layer  stiffness 
properties  E,,  E2,u12,  G12  and  the  layer  geometry  and  orientation. 
Similarly,  can  be  determined  from  the  damping  characteristics  of 
each  layer  and  the  partiuclar  viscoelastic  model  that  is  adopted. 

There  are  numerous  papers  that  develop  the  laminae  stiffness 
properties  in  terms  of  the  micromechanical  properties  Em,  Ef,  vml  u, 
and  the  volume  fraction  v  of,  say,  fibers.  Virturally  all  studies 
adopt  the  law  of  mixtures [1]  to  obtain  the  principal  Young's 
modulus  E,  and  Poisson's  ratio  vn  in  terms  of  the  aforementioned 
micromechanical  properties.  There  are  some  other  methods  that 
employ  elasticity  and  its  theorems  on  minimum  potential  and 
complementary  energy  to  obtain  bounds  on  E,  and  uI2  [2-4]. 
Fortunately,  these  upper  and  lower  bounds  are  fairly  close  to  each 
other.  We  have  found  by  direct  calculation  that  the  frequency 
response  of  angle-ply  laminates  is  relatively  insensitive  to  the 
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various  assumed  values  of  E,  and  uI2  between  their  upper  and  lower 
bounds . 

The  stiffnesses  E2  and  G12  can  vary  greatly  depending  upon  the 
assumed  modelling  assumptions.  Likewise,  elasticity  approaches 
show  that  there  are  wide  percentage  variations  between  their 
permissible  maximum  and  minimum  values.  Howerver,  for  most  fiber- 
reinforced  laminae,  these  stiffnesses  are  small  compared  to  E,.  As 
a  result,  the  frequency  response  of  angle-ply  laminates  is 
relatively  independent  of  the  values  E2  and  GI2,  provided  these 
values  are  selected  within  their  elasticity  defined  limits. 
Consequently,  insofar  as  the  objective  to  determine  the  frequency 
response  of  laminated  composites  in  terms  of  its  micromechanical 
elastic  properties  is  concerned,  it  appears  that  the  result  is 
relatively  independent  of  the  micromechanical  model. 

Viscoelastic  modelling,  however,  is  a  different  matter.  In 
contrast  to  its  elastic  counterpart,  there  are  very  few  published 
papers  on  this  subject.  In  one  such  paper  [5],  some  viscoelastic 
properties  are  dtermined  both  experimentally  and  analytically. 

A  second  major  objective  of  this  research  was  to  develop  new 
models  which  retain  the  local  discontinuities  of  multiphase 
materials  and  are  capable  of  characterizing  composite  fiber- 
reinforced  laminates.  The  approach  was  to  adapt  a  methodology, 
called  the  load  correction  method  (LCM)  [6]  for  the  static 
analysis  of  repetitive  lattice  structures  having  local 
discontinuities  [7]  to  the  dynamic  analysis  of  composites. 
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Chapters  2-5  are  deveoted  to  explicity  obtaining  the  solution 
for  the  frequency  of  laminated  symmetric  plates  in  terms  of  the 
laminates  stiffness  and  damping  properties.  Chapter's  6  and  7 
report  on  the  progress  i".  developing  the  LCM  approach. 

Specifically,  Chapter  2  is  based  upon  a  paper  in  which  the 
distributed  parameter  form  of  the  eigenvalue  equation  for  laminates 
and  sensitivity  analysis  is  used  to  develop  an  explicity  solution 
for  the  frequency  response  of  elastic  laminates  subjected  to 
various  boundary  conditions.  This  paper  was  presented  at  OPTI  '89, 
an  international  conference  held  in  Southampton,  UK. 

In  Chapter  3,  the  general  sensitivity  analysis  method 
presented  by  Reiss  [7],  is  specialized  to  discrete  self-adjoint 
eigenvalue  problems.  This  material  was  presented  in  Valparaiso, 
Chile,  1991  at  the  2nd  Pan  American  Congress  on  Applied  Mechanics. 

Chapter  4  extends  the  Pi's  earlier  results  to  include  viscous 
damping.  Sensitivity  analysis  is  again  used  for  both  distributed 
and  discrete  models.  Specific  accurate  approximations  for  the 
complex  eigenvalues  are  given  in  terms  of  the  stiffness,  mass  and 
damping  matrices.  This  chapter  was  presented  at  the  ASME  Design 
Automation  Conference  in  Chicago,  Septermber,  1990. 

Chapter  5  employs  the  methods  of  Chapter  3  to  obtain  still  new 
solutions  for  the  frequency  response  of  laminated  plates.  This 
solution  is  slightly  more  accurate  than  the  one  obtained  in  Chapter 
2.  Moreover,  unlike  the  method  in  Chapter  2,  this  approach  is 
rigorous  and  can  be  used  to  quantify  the  error  in  the 
approximation.  These  results,  originally  intended  for  presentation 
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at  the  cancelled  IMAC  meeting,  were  presented  by  invitation  at  the 
Florence  Modal  Analysis  Conference,  Florence,  Italy,  in  September, 
1991. 

Chapter  6  summarizes  the  development  of  the  load  correction 
method  to  date.  Some  of  the  computational  aspects  of  this  method 
are  presented  a  Chapter  7,  a  paper  also  presented  at  PACAM  II. 

Chapter  8  summarizes  some  of  the  conclusions  of  this  report 
and  makes  some  specific  recommendations  for  further  study. 
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2.  The  Fundamental  Frequency  of  Symmetric 

Laminates  Determined  by  Gigensensitivity  Analysis 
Bo  Qian  and  Robert  Reiss 

Department  of  Mechanical  Engineering,  Howard  University, 
Washington,  D.C.  20059,  USA 

ABSTRACT 

A  new  closed-form  approximate  solution  for  the  fundamental 
frequency  of  symmetric  rectangular  laminates  is  presented  for 
all  classical  combinations  of  clamped  and/or  hinged  edge  sup¬ 
port  conditions-  The  derivation  consists  of  an  expansion  of 
the  fundamental  frequency  in  n  Mucl.iurin  series  in  the  cou¬ 
pling  stiffnesses  Difc  and  035  and  then  truncating  after  qua¬ 
dratic  terms.  The  r<  .tired  sensitivity  derivatives  are  calcu¬ 
lated  using  a  method  developed  by  Reiss  [lj.  Calculated  re¬ 
sults  obtained  from  the  new  approximate  formula  are  generally 
within  27.  of  the  frequency  obtained  from  conventional  numeri¬ 
cal  methods,  but  require  but  a  fraction  of  their  computational 
effort. 

INTRODUCTION 

An  accurate  method  to  determine  the  flexural  fundamental  fre¬ 
quency  of  symmetric  rectangular  laminates  is  necessary  for  mod¬ 
ern  laminated  composite  ^°sign.  Closed-form  exact  solutions 
are  available  only  for  simply-supported  orthotropic  laminates 
[2].  Closed-form  approximate,  but  very  accurate,  solutions  for 
orthotropic  plates,  none  of  whose  edges  are  free,  have  been 
obtained  by  Hearmon  [3]  using  a  one-term  Ray le igh-Ritz  proce¬ 
dure  . 

The  literature  c  ■ jo  contains  a  few  closed-form  approxima¬ 
tions  for  anisotropic  laminates.  Bert  considered  general  sym¬ 
metric  rectangular  laminates  whose  edges  are  either  simply- 
supported  [4]  or  clamped  [5].  He  assumed  a  linear  relation¬ 
ship  between  eigenvalue  and  flexural  stiffnesses,  and  deter¬ 
mined  the  unknown  constants  from  numerical  calculations.  A 
different  approach  was  taken  by  Reiss  and  co-workers.  Reiss, 
Ramachandran  and  Qian  [6],  using  exact  sensitivity  derivatives, 
obtained  an  approximate  solution  for  the  natural  frequencies 
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of  a  four-layer  symmetric  simply-supported  angle-ply  laminate. 

The  frequencies  . . I  were  |»  re  setii  etl  in  r  he  form  of  a  trun- 

C’.ed  Fourier  series  in  the  ply  angle.  Reiss  and  Qian  [7] 
chen  generalized  those  results  to  include  all  simply-supported 
symmetric  laminates. 

The  purpose  of  this  paper  is  to  develop  one  simple  formu¬ 
la  to  determine,  app rox ima t e 1 y ,  the  fundamental  frequency  of 
symmetric  lammaces  subject  to  all  combinations  of  hinged  and 
clamped  edge  conditions. 

PROBLEM  STATEMENT 


The  eigenvalue  equation  for  a  freely  vibrating  symmetric  lami¬ 
nate  is 


11 


w, 


xxxx 


D,«  w, 

22  yyyy 


2<D12  . 


2D.,)  "• 

6o  xxyy 


4  D,,  w, 

16  xxxy 


4  D 


26 


w, 


xyyy 


=  0  hi 


(1) 


where  w  denotes  the  mode  shape,  w  is  the  corresponding  fre¬ 
quency,  p  is  the  mass  density,  are  the  flexural  stiff¬ 

nesses  [l],  x  and  y  are  the  coordinates  in  the  plane  of  the 
laminate,  and  indices  following  a  comma  denote  differentiation 
with  respect  to  the  indicated  argument. 


•  Each  of  the  edges  x=0  and  x=a  is  either  hinged  or  clamped. 
The  appropriate  boundary  rnndii  tone  are 


hinged : 


l  1 


2  It 


10 


xv 


=  0. 


clamped : 


w  =  w  ,  O 

X 


\2) 


Similarly,  each  edc,'>  v~0  and  y--h  is  also  either  hinged  or 
clamped.  Titus  these  boundary  conditions  are 

hinged:  w  =  2D26  *  D22  w,yy  =  0 

(3) 

clamped:  w  =  w,^  =  0 

The  objective  of  this  paper  is  to  determine  the  lowest 
frequency  for  each  of  the  combinations  of  boundary  condi¬ 
tions  (2)  and  (3)  in  terms  of  the  material  parameters  Djj  and 
plate  aspect  ratio  R: 

R  *  a  /  b  (  4 ) 

EIGENVALUE  DIFFERENTIATION 

Before  proceeding  with  the  development  of  the  solution  to  the 
stated  problem,  it  is  useful  to  review  some  general  mathemati- 
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cal  formulae  that  are  essential  to  the  current  approach. 

Consider  the  class  of  eigenvalue  problems  which  can  be 
represented  in  the  form 

T*£(S)T  w  =  \  M  w  (5) 

mn 

to  which  appropriate  mixed  boundary  conditions  must  be  append¬ 
ed.  A  detailed  discussion  of  the  requirements  the  operators 
appearing  in  Equation  (5)  and  the  boundary  operators  must  sat¬ 
isfy  can  be  found  in  the  paper  by  Reiss  and  Haug  [8].  For  our 
present  purposes,  it  suffices  to  note  that  T  and  T*  are  L2 
adjoint  differential  operators,  S  is  a  collection  of  variable 
parameters,  and  E  and  M  are  positive  stiffness  and  mass  opera¬ 
tors,  respectively.  The  double  indices  m  and  n  are  included 
in  order  to  facilitate  a  solution  for  the  two-dimens ional  lam¬ 
inate. 

The  eigenvalue  Xmn  is  a  functional  of  the  parameters  de¬ 
noted  symbolically  by  S.  And  if  these  parameters  change  by  a 
small  amount  6S,  Xmn  changes  by  a  small  amount  6Xmn  [ 1  ] : 

6Xmn  =  (T  w  ,  6  E  T  w  )  (6) 


where  the  normalized  eigenfunctions  satisfy 

(w  ,  M  w  )  =  1  (7) 

In  Equations  (6)  and  (7),  (•,  •)  denotes  the  usual  L2  inner 
product  over  the  domain  spanned  by  the  functions  involved  -  in 
the  case,  the  planar  area  of  the  plate.  Furthermore,  the 
second  variation  of  the  eigenvalue  is  given  by  the  expression 
[1]« 

S2X  =  (T  wmn,  «2E  T  wmn) 
mn 

-2  I  l  wlj,6  E  Twmn)2  (8) 

i,j  ^m,n  *ij  -  *mn 

In  Equation  (9),  the  summation  takes  place  over  all  pairs  of 
positive  integers  i,j  except  i=m  and  j=n.  Equations  (6)  and 
(8)  are  valid  provided  the  eigenvalue  xmn  is  not  repeated. 

Similarly,  the  first  variation  of  an  eigenmode  wmn  asso¬ 
ciated  with  a  non-repcated  eigenvalue  1  ■, 


.  mn  ,  „  w  J 

«w  =  -  l  l  - - — 

i  ,  j  ^m ,  n  i  j  -  mn 


(T  wlj,6  E  T  wmn)  (9) 
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1c  is  clear  from  Equations  (8)  and  (9),  that  it  is  necessary 
to  know  the  complete  set  of  eigenvalues  and  eigenfunctions  for 
s  erne  pa  ram-  ;  e  :  :•  in  !  .  >  c  a  !  .  1 1  I  a  f  c  I  hr  r.ni.il  ii'im  ■  \ 

.  mil 

and  Oi  w‘ml  exactly. 

For  symmetric  laminated  plates,  we  can  identify 


T  =  [3  /3x  3  /3y  2  32/3x3y] 


(10) 


T*  =  T  ,  M 


X  =  a) 
mn  mn 


For  the  specific  boundary  conditions  (2)  and  (3),  we  can  se¬ 
lect  the  operator  E  as: 


(ID 


where 


'Dll 

D 

< 

E  = 

D 

< 

°22 

°26 

-D16 

°26 

0  _ 

D  = 
< 

°12  +  2 

°66 

(12) 


Therefore  the  eigenvalues  are  functions  of  the  five  material 
parameters  Dn,  D22.  D* ,  and‘D26. 

TH*:  ORTHOTROPIC  SOLUTION 


For  orthotropi.  I.nnin.ii  i”.,  i  In-  coupling  stiffnesses  D15  and 
vanish.  A  lompli-ii-  •.«•(  nl  «• i  gcnv.i  I  ues  and  eigenfunctions  can 
be  obtained  1 1  1)<  also  vanishes.  Thus  wc  identify  S  with  D< 

and  consider 


U 


11 


~llltl 

W, 

xxxx 


22 


miiii 

w’yyyy 


p 


x 

mn 


~mn 

w 


(13) 


subject  to  tho  app  rop  r  i  a  r  1-  boundary  conditions  (2)  and  (3)  as 
simplified  lor  an  orihotropic  plate.  The  solution,  obtained 
by  separating  variables,  is 


aIHII 

w 


X  (-)  Y 
m  a  n 


(14) 


P 


X 

mn 


D„  v 
22  n 


(15) 


where  Xm  and  Yn  are  the  normalized  beam  shape  functions  [3]  and 
and  their  corresponding  frequencies,  i.e. 


(16) 
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IV 

In  Equation  (16),  (  )  denotes  the  fourth  derivative  with  re¬ 
spect  to  the  argument  of  the  indicated  function. 

The  solution  for  the  orthotropic  laminate  can  be  express¬ 
ed  as  a  Maclaurin  series  in  D  : 

ic 


a£ 


inn 


=  £ 


mn 


mn 


3D. 


Ok  *  k 


32  £ 


3  D  - 


mn  „  2 
2  D* 


(17) 


where  (''l  denotes  evaluation  of  the  indicated  function  at 
Dk=0.  The  linear  term  is  calculated  by  observing  that 

'0  1  O' 


6E  = 


3D. 


0 

0 


(18) 


and  then  substituting  this  result  together  with  Equation  (14) 
into  Equation  (6)  to  obtain 


a£mn/3D, 


2  2 

=  2  C  c  /p  a  b 

mm  nn 


where 


Cij  =-  (Xi’  Xi) 


c.  .  =  -  (Y!  ,  Y' ) 
AJ  1  J 


(19) 

(20) 


The  second  derivative  is  also  readily  obtained  by  observing 

that  62E  =  a2E/3DK2  vanishes  identically.  Now,  substitution 

of  Equations  (14),  (15)  and  (18)  into  Equation  (8)  vields 

2 ~  C.'  2  c  .2 

,  3  Xmn  ,  _  _ _ lm  nj  _ 

^  1  ~  “  i,j  *m,n 


3  D. 


,4  4.  4  .  4  4  4 

D  (u.  -p  )b  +D._ ( v  .  -v  )a 
111  m  22j  n 


(21) 


The  orthotropic  natural  frequencies  Xmn  become 


0  A 


mn 


11 


4  4 

%  'a 
m 


22 


4.,  4 
v  /l) 
n 


2  C 


V  • 


-  4 


_  2  2  _  2 
C .  c  D. 
lm  n  j  * 


^  4  4  b  -r 

i,j  ^m,n  D. , (p .  -p  )b  +D,_(v.  -v  ’)a 
Him  Z2j  n 


— 4  *°<V> 


(22) 


Equation  (22)  is  exact  to  within  terms  of  order  DK  .  The 
first  three  terms  on  the  right-hand  side  of  Equation  (22)  is 
precisely  the  solution  obtained  by  Hearmon  [3]  using  an  entire¬ 
ly  different  approach.  The  quadratic  term  in  Equation  (22) 
represents  a  correction  to  Hearmon's  solution.  Routine  calcu¬ 
lations  for  a  wide  range  of  parameters  .  and  plate  aspect 
ratios  R  show  that  the  second  order  correction  is  not  only 
negligible  for  the  fundamental  frequency  of  orthotropic  plates, 
but  also  for  higher  natural  frequencies.  The  only  exceptions 
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occur  for  higher  frequencies  and  select  values  of  R,  namely 
those  for  which  the  denominator  of  the  quadratic  terms  be¬ 
comes  small  for  some  pair  of  indices  i  and  j.  This  corre¬ 
sponds  to  Lite  when  I'qn.i!  inn  (13)  admits  repeated  eigen¬ 

values,  that  is,  when  the  right-hand  side  of  Equation  (15) 
has  the  same  value  for  two  distinct  pairs  (m,n). 


Using  the  foregoing  approach,  the  eigenfunctions  w  for 
the  orthotropic  laminate  can  be  shown  to  satisfy 

wmn  =  X  Y  -D  .  E.  ,£  a..  X  .  Y  .  +  . . .  (23) 

m  n  k  i,j  *m,n  ijmn  l  j 

where 


2  C.  c  . 
lm  n  j 

U  4  -2  ”  "  4  4  2 

D. , (m  .  -u  )R  +  D0_ ( v  .  -v  )R 
11  l  m  22  j  n 


A  few  observations  are  in  order.  For  most  cases  that  we  cal¬ 
culated,  the  quantity  aijmn  D<  was  small  compared  to  unity 
for  all  i  and  j.  In  such  cases,  the  solution  w*1111  is  a  suffi¬ 
ciently  accurate  approximation  to  w"11*.  Furthermore,  in  cases 
when  wnm  cannot  approximate  wran,  the  terms  ajjmn  D„.  was  found 
to  be  comparable  to  unity  only  for  one  set  of  indices  i  and  j; 
more  importantly,  this  term  contributed  a  negligible  amount  to 
the  inner  products  required  for  subsequent  calculations.  Thus, 
throughout  the  remainder  of  this  investigation,  we  will  sub¬ 
stitute  w"10  for  wmn  in  all  inner  products. 

THE  ANISOTROPIC  SOLUTION 


Like  the  orthotropic  solution,  the  anisotropic  solution  is 
developed  by  expressing  it  in  a  Maclaurin  series.  Here,  the 
parameters  S  are  identified  with  Di6  and  D26.  Thus 


where  the  symbol  (  )denotes  evaluation  of  the  indicated  func¬ 
tion  at  the  orthotropic  solution  Di6=D26=0. 

Since  the  orthotropic  solution  is  now  completely  known, 
albeit  approximately,  each  of  the  coefficients  in  the  series 
(25)  can  be  evaluated.  These  calculations,  when  substituted 
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back  inco  Equation  (25),  lead  to  the  desired  closed-form  ap¬ 
proximation  to  the  fundamental  frequency,  namely 

»bSll  ■  °U  *  D22  “l4  *  2  D<  CHc11/b2 


,-2 


2.  ,2 


4  R-  Xl[(Hil8lj+Hligjl)D16  +  (G1.h.rG.1hlj)R^D26] 


^i^V>-°22<vv>4* 


(26) 


where 


2D.  82(CllcJJ  -  Cucn)r‘ 


Gu  - (x;-  v 

«u  -  <*;• 

v.) 

j 

Hu  ■  <xi  •  x;> 

hij  ■  <*•;• 

Y.) 

J 

(27) 


Apart  from  the  orthotropic  approximations  employed  in  this 
derivation,  Equation  (26)  is  correct  to  within  quartic  terms 
in  the  coupling  stiffnesses. 

Inclusion  of  higher  order  terms  in  the  Maclaurin  series 
is  unlikely  to  yield  any  significant  improvement  in  the  accu¬ 
racy  of  the  approximation  (26)  unless  the  orthotropic  solution 
is  also  improved.  The  only  exception  is  the  simply-supported 
plate  for  which  the  beam  shape  functions  provide  the  exact 
orthotropic  solution.  Nevertheless,  even  for  this  plate,  the 
addition  of  higher  order  terms  negates  the  simplicity  of  the 
approximation.  Furthermore,  as  the  numerical  illustrations  in 
Tables  I-I1I  demonstrate,  Equation  (26)  is  remarkably  accurate. 

NUMERICAL  EXAMPLES 

In  the  following  tables,  the  validity  of  Equation  (26)  is 
established  by  comparing  the  fundamental  frequency  calculated 
from  it  with  the  fundamental  frequency  obtained  numerically. 
Since  Equation  (26)  uses  the  vibrating  beam  shape  functions, 
the  numerical  approach  selected  was  the  Ritz  method  using  these 
same  shape  functions.  Thus  the  fundamental  mode  is  expressed 

11  N  N 

w  =  ih  j=i  °ij  Xi(x/a)  Yj(y/b>  (28) 

It  is  convenient  to  introduce  the  following  non-dimension¬ 
al  fundamental  frequency. 
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k  =  Wjj  •'Wu  j 

where  is  an  invariant  material  property  [2]  and  h  is  the 
common  thickness  of  the  plies. 

Results  are  presented  for  balanced  symmetric  four-ply  plates 
[ 9/— ©/— 0/ 6]  -  Each  table  includes  results  for  typical  boron 
epoxy  (B/E)  and  high  modulus  graphite  epoxy  (G/E)  defined  by 
the  following  data. 

B/E  G/E 


E 1/E2  10  40 

Ex/G12  40  80 


\> 


12 


.30 


.25 


The  tables  have  been  developed  for  different  boundary 
conditions,  whose  nomenclature  is  S  for  a  simply  supported 
edge  and  C  for  a  clamped  edge.  The  first  letter  denotes  the 
support  condition  on  the  edge  x=0  and  the  subsequent  letters 
are  for  edges  proceeding  counterclockwise  around  the  plate. 

For  example  CCCS  means  clamped  on  all  sides  except  y=b  where 
it  is  simply  supported.  The  columns  labelled  tci  are  calcula¬ 
ted  by  the  Ritz  method  with  N=ll  i.e.  121  terms  in  Equation 
(28).  In  most  cases,  the  Ritz  method  has  converged  with  N=ll; 
in  those  cases  for  which  convergence  is  not  evident  for  N=ll, 
comparison  with  the  N=7  and  N=9  term  solutions  suggests  that 
the  indicated  value  is  within  1%  of  the  correct  frequency. 

It  should  be  observed  that  in  all  cases  presented  the 
approximate  solution  (26)  is  within  27.  of  the  Rayleigh-Ritz 
solution,  and  often  much  better  than  that.  Calculations  for 
other  boundary  conditions  and  other  material  data  do  not  al¬ 
ter  this  conclusion. 
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TABLE  I: 

RESULTS 

FOR  CCCC  PLATE 

9° 

B/E 

kl 

G/E 

kl 

kl 

0 

84.17 

84.19 

R=1 

84.28 

84.30 

15 

82.26 

82.38 

81.59 

81.79 

30 

78.65 

79.19 

76.62 

7  7.49 

45 

77.00 

77.76 

74.46 

75.59 

0 

25.28 

25.28 

R=5 

13.64 

13.65 

30 

31.51 

31.61 

25.74 

25.86 

60 

61.21 

61.26 

62.95 

62.99 

90 

78.65 

78.66 

82.54 

82.54 

TABLE 

II j  RESULTS 

FOR  SCSC  PLATE 

B/E' 

G/E 

9° 

kl 

kl 

*1 

kl 

R=1 

0 

80.82 

80.83 

83.33 

83.33 

30 

74.27 

74.40 

75.47 

75.62 

60 

57.36 

57.96 

55.05 

55.93 

90 

45.33 

45.36 

39.96 

39.97 

R=5 

0 

11.81 

11.81 

6.94 

6.94 

30 

15.67 

15.87 

13.44 

13.65 

60 

27.64 

27.73 

28.33 

28.42 

90 

34.80 

34.81 

36.46 

36.46 
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TABLE  III:  RESULTS  FOR  CCCS  PLATE 


B/E 

C/E 

0° 

*1 

kl 

kl 

kl 

R=1 

0 

61.84 

61.87 

59.34 

59.36 

30 

64.70 

65.44 

61.54 

62.69 

60 

75.01 

75.42 

74.23 

75.27 

90 

82.10 

82.13 

83.70 

83.72 

R=5 

0 

25.16 

25.16 

13.42 

13.42 

30 

31.19 

31.50 

25.64 

24.76 

60 

61.19 

61.24 

62.94 

62.99 

90 

78.64 

78.65 

82.53 

82.54 
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ABSTRACT 

A  new  method  is  presented  to  determine  approximate  closed- 
form  solutions  for  the  eigenvalues  of  self-adjoint  structures  with 
non-uniform  stiffness  and  mass  distributions.  The  distributed 
parameter  eigenvalue  equation  is  cast  into  an  equivalent  discrete 
formulation  with  an  infinite  number  of  degrees  of  freedom.  The 
resulting  stiffness  and  mass  matrices  are  each  decomposed  into  two 
matrices,  one  diagonal  and  the  other  with  zero  diagonal  elements. 
Sensitivity  derivatives  are  used  to  expand  the  eigenvalues  in  a 
power  series  of  the  zero  diagonal  matrices. 


INTRODUCTION 

Eigensensitivity  methods  have  been  ased  primarily  in  optimum 
design  or  reanalysis  of  structures.  Recently,  these  methods  have 
been  used  to  obtain  approximate  closed-form  solutions  to  specific 
eigenvalue  problems  [1] .  In  this  paper,  eigensensitivity  analysis 
is  used  to  obtain  a  simple  expression  which  approximates  the 
eigenvalues  for  a  wide-class  of  self-adjoint  eigenvalue 
problems . 


Consider  the  eigenvalue  problem  in  abstract  form 

T*  E(S)  T  Ui  -  M <S)  uj  (1) 

subject  to  appropriate  mixed  boundary  conditions.  Here,  T  and  T’ 
are  L  2  adjoint  differential  operators,  E  is  the  positive,  symmetric 
stiffness  operator,  M  is  another  positive  symmetric  operator  and 
S  denotes  a  collection  of  parameters  upon  which  both  E  and  M 
depend.  Furthermore  is  the  eigenvalue  with  associated 

eigenfunction  u^. 
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There  are  many  ways  to  represent  Eq.  (1)  by  an  equivalent 
discrete  system.  In  this  paper,  the  Rayleigh-Ritz  method  will  be 
used.  If  {  is  a  complete  kinematically  admissible  set  of  basis 

functions,  and 

u.  =  I  a  .  lit  (2) 

x  n  m  n 

then  stationarity  of  the  Rayleigh  quotient  requires 

ko.  -  x^a.  O) 

where  the  components  Kjj  and  Mj.  of  the  stiffness  matrix  K  and  the 
matrix  H  respectively  satisfy 

-  (T<|>.,  ETip.)  (4) 

M.  .  =  (  <i>  - ,  .  ) 

ij  i  J 

Here,  (  •  ,  •  )  denotes  the  L^- inner  product.  Also  ain  is  the  n-th 
component  of  the  infinitely  dimensional  eigenvector  q  .  . 


FORMAL  SOLUTION 

It  is  desired  to  solve  Eq.  (3)  for  the  eigenvalues  Xj.  Toward 
this  end  it  is  convenient  to  introduce  parameters  S  and  S2  and  the 
one  parameter  family  of  matrices  K  (S^)  and  M  (S2)  aefined  by 

K  (SL)  -  Kd  +  SjAr  (5) 

M  (S2  )  -  Md  +  SjflM 

In  Eq.  (5),  Kd  and  MD  are  diagonal  matrices  whose  components 
consist,  respectively,  of  the  diagonal  elements  of  K  and  M;  also 
AK  and  AM  have  zero  diagonal  elements,  and  off  diagonal  components 
equal  to  those  of  K  and  M,  respectively.  Evidently, 


K  (1)  -  K 


M  (1)  =  M 


Consider  the  eigenvalue  equation 

ai  (Srs2)  «  x.<s1,s2)  k  " 1  ( s  p  M-1  (S  2)  at  <s  ,s2)  p) 

A  -  1  -  1 

If  S2*  0,  K  and  M  are  diagonal,  and  therefore  Eq.  (7) 

readily  admits  the  solution 


X.<0,0)  =  Ku  MU 


a. (0,0)  =  m:j  e 

i  hi 
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where  a.  is  the  unit  Cartesian  base  vector. 

i 


Equation  (7)  is  a  special  case  of  the  eigenvalue  equation 
discussed  by  Reiss  (2]  in  which  an  explicit  representation  for  the 
variations  of  the  eigenvalues  and  eigenvectors  was  determined.  In 
the  context  of  Eq.  (7.',  these  variations  become 

<5X  .  =  a.  AKa.  6  s.  -  x  a.  am  a.  i  s, 

i  i  i  i  ii  i  - 


6 


2  I 

nH 


(xianAMai6s2-anAicai  6  s^2 


2  6  X 


a  a  Ma.  6  s 

i  i 


provided  the  eigenvectors  <3i  are  normalized  with 
Furthermore,  the  variations  (9)  and  (10)  satisfy 


(9) 


(10) 


A 

respect  to  M. 


\  i(S1+6S1,S2+  6S2)  -  X  .(S^)  =  6X  .  +  |  6“X  i  +  ...  (11) 

Since  (1,1)  =  X4  ,  the  solution  for  may  be  obtained  by 

substituting  Eqs.  (9)  and  (10)  into  Eq.  (11)  and  evaluating  at 
S  j  ■  S2  m  0,  6S2  »  <5  S  2  =  1.  Thus  if  only  second  order  terms  in 
6 Sj  are  retained. 


„  (M.  .AK  .  -  K.  AM  .Y 

-1  -2  _  11  ni  11  m 

X  .  »  K..M.!  -  M  l  - 

1  ii  ii  11  n^i  m.  .  k  -  M  K. 


11  nn 


nn  11 


(12) 


Even  higher  order  variations  could  be  calculated,  however,  only  at 
the  expense  of  the  simplicity  of  the  representation  (12) . 


EXAMPLE 

Consider  the  torsional  vibrations  of  a  fixed-free  shaft  whose 
radius  r  varies  with  the  axial  dimensions  x.  The  non-dimensional 
eigenvalue  equation  is 

A  if  \  4 

-  (r  (x)  u  )  =  A  r  (x)u 

l  '  ... 

where  u(0)  *  u  (1)  =  0  and  (  )  =  d/dx.  A  convenient  basis  is 

^n*sin(n~2)TTX  ,  n=l,  2,  .  .  . . 

An  approximate  solution  for  X  ^an  be  obtained  from  Eqs.  (4)  and 
(12)  for  any  r(x).  Assuming  r“=  a-x,  a  >  1,  the  approximate 
eigenvalues  obtained  from  Eq.  (12)  are  tabulated  in  Table  1  for 
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two  different  values  of  a.  Note  that  as  a-*®  ,  the  shaft  tends  to 
a  uniform  one  and  4>n  approaches  the  exact  eigenfunctions.  The 
first  column  is  the  zeroth  order  solution  Kii  obtained  by 

ignoring  the  off-diagonal  elements  of  M  and  K.  The  next  column 
represents  the  second  order  correction  obtained  by  adding  i<5 
to  the  first  column.  Finally,  since  this  tapered  shaft  admits  an 
exact  solution,  the  last  column  contains  the  exact  eigenvalues. 


Table  1.  Approximate  eigenvalues  for  shaft 


i 

a  *  1 . 5 

a  =  3 

1 

5.417 

5.211 

5.239 

3.40/ 

3.371 

3.373 

2 

24.13 

25.24 

25.88 

23.01 

23.17 

23.19 

3 

63.56 

63.96 

65.55 

62.48 

62.60 

62.68 

4 

122.8 

123.1 

124.8 

121.7 

121.9 

121.9 

5 

201.7 

202.1 

203.8 

200.7 

200.8 

200.9 

Other  illustrations,  such  as  the  vibration  of  non-uniform  beams 
and  buckling  of  bars  with  variable  cross-section  can  be  solved  with 
similar  ease  and  accuracy. 
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ABSTRACT 

A  new  method  is  presented  to  determine 
approximate  closed-form  solutions  for  the 
c omplex- valued  frequencies  of  moderately  damped 
linearly  elastic  structures.  The  approach  is 
equally  applicable  to  finite  degree  of  freedom 
systems  and  distributed  parameter  systems.  The 
damping  operator  is  split  into  two  components, 
the  first  of  which  uncouples  the  quadratic 
eigenvalue  equation,  and  the  eigenvalues  are 
expressed  as  a  power  aeries  in  the  second 
component  of  the  damping  operator.  Specific 
numerical  examples  include  both  finite  degree 
of  freedom  and  distributed  parameter  systems. 

It  is  shown  that  for  moderate  damping,  that  is, 
when  the  second  component  of  the  damping 
operator  is  small,  but  not  negligible,  the 
series  solution  truncated  after  quadratic  terms 
provides  an  excellent  approximation  to  the  true 
eigenvalues . 


Kuseyin  (1973),  are  among  the  first  to  develop 
explicit  closed-form  expressions  for  the  design 
derivatives  of  self-adjoint  and  non  self- 
adjoint  systems,  respectively.  Sensitivity 
derivatives  for  distributed  parameter  systems 
were  presented  recently  by  Reiss  (1986)  by 
casting  the  governing  self-adjoint  differential 
equation  into  integral  form  using  Green's 
function. 

In  this  study,  the  sensitivity  derivatives 
of  the  complex  eigenvalues  associated  with 
viscously  damped  structures  are  explicitly 
determined.  Explicit  dependence  of  the 
eigenvalues  and  eigenfunctions  upon  boundary 
conditions  is  eliminated  by  presenting  the 
quadratic  eigenvalue  equation  in  integral  form 
For  simplicity,  it  is  assumed  that  the  design 
variables  affect  the  eigenvalues  and  eigen¬ 
functions  only  through  the  damping  operator. 

PROBLEM  STATEMENT 


INTRODUCTION 

Eigensensitivity  analysis  is  an  exception¬ 
ally  useful  tool  for  the  reanalysis  of  eigen¬ 
values  and  eigenvectors.  Although  it  is 
essentially  a  perturbation  technique,  such 
methods  have  been  used  successfully  to 
determine  optimum  eigenvalue  designs  as  well  as 
approximate  closed-form  expressions  for  the 
eigenvalues  of  structural  systems. 

In  eigensensitivity  techniques,  the 
eigenvalues  and  eigenvector  are  assumed  to  be 
piecewise  differentiable  functions  of  some 
specified  set  of  parameters,  called  design 
parameters.  It  is  the  derivatives  of  the 
eigenvalues  and  eigenvectors  with  respect  to 
the  design  parameters  that  are  sought.  Tor 
systems  with  a  finite  number  of  degrees  of 
freedom.  Fox  and  Kapoor  (1968)  and  Plaut  and 


The  closed,  bounded  and  regular  domain  of 
an  elastic  structure  is  denoted  by  n  ♦  ,.i  where 
the  interior  of  the  domain  is  ft  and  its 
boundary  is  3ft.  Furthermore,  3!)  is  the  union 
of  two  regular  complementary  subsets  -  and 
Jft2 .  A  typical  point  in  ft  +  3ft  is  denoted  by  x 
and  a  typical  time  is  denoted  by  t.  The  design 
variables,  which  may  depend  upon  x  but  not  t, 
are  collectively  denoted  by  Six)  . 

The  unforced  motion  of  a  linear  elastic 
structure  vibrating  in  a  viscous  medium  (dis¬ 
tributed  damping)  can  be  represented  by  the 
equation 

T*ETu  +  C  (S)  u  ♦  Mu  -  0  (x,t)  tftx  (0,-)  (la) 

subject  to  the  mixed  boundary  conditions 

B  u  -  0  (X,  t>  e  3(3;  X  (0,  -  )  (lb) 
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B*ET  u  -  0 


<x,  t)  c 


>02*  JO-  -> 


(lc)  DESIGN  DERIVATIVES 


Hare,  T  and  T*  are  Lj  adjoint  differential 
operators,  B  and  B*  are  the  corresponding 
adjoint  boundary  operators,  and  the  positive, 
symmetric  operators  E,  M  and  C  are,  respec¬ 
tively,  the  stiffness,  mass,  and  damping 
operators.  The  response  u  depends  upon  both  x 
and  t,  and  the  dot  above  a  symbol  indicates  its 
time  derivative. 

The  eigenvalue  problem  under  consideration 
is  obtained  by  substituting  the  separable  form 

u(x,t)  -  U  (x)  eJt 

into  Eqs .  (la,b,c).  This  results  in  the 
quadratic  eigenvalue  equation 

T*ET  Ui  -  -«i  C(S)Oi  -  U.J2  M  Ui  x  c  n  (2a) 


The  eigenvalue  equation  (3)  is  not  self- 
adjoint  and,  therefore,  the  eigenfunctions  are 
not  orthogonal  in  any  useful  sense.  Following 
Plaut  and  Huseyin  (1973),  Eq.  (3)  car.  be  re¬ 
duced  to  standard  form  by  introducing  v.  -  .  u 
and  the  identity  operator  6  with  respect  to 
(  •  ,  •  )fi  ■  Thus  Eq.  (3)  becomes 

«i  -  (G,  C  Ui  +  M  VJjj 

(6) 

V,  ■  u  ■  (  4  U  \ 
i  i  '  '  ui'n 

In  matrix  form,  Eq.  (6)  is 


0 

■r4 

3 

1 

1 

c 

°) 

(c 

M) 

C") 

'Vi' 

Vo 

<5«  ' 

'  - 1 

o  ' 

W/  (7) 

8  Ux  - 

0 

x  c  jflx 

(2b) 

B*ET  Ux  - 

0 

x  c  3(!2 

(2c) 

In  Eq.  (2),  w1  is  the  complex-valued  eigenvalue 
and  Ui  is  the  associated  eigenfunction.  Fol¬ 
lowing  Reiss  (1986),  Eqs.  (2)  are  recast  in  the 
simpler  single  integral  equation 

U.  (y)  -  -wj.  (G  ( • ,  y) ,  C(S)  OJn"  ui^G( -,y) ,  MUi)n 

(3) 


where  the  operators  G-  and  are  defined 
through 


G- 

5 

(G  , 

•  ) 

4  • 

3 

(  5  , 

•  ) 

In  Eqs.  (8),  G  is  the  same  self-adjoint  Green's 
function  appearing  in  Eq.  (3) . 

The  first  variation  of  Eq.  (7)  is  now 
desired.  Thus  it  is  assumed  that  S  changes  by 
4S,  so  that  the  first  variation  of  Eq. (7)  is 


where 

(a,b)n  •'  !a  a  b  dx 

and  G  is  the  Green's  function  associated  with 
the  operators  on  the  left  side  of  Eqs.  (2),  and 
a  is  the  complex  conjugate  of  a. 


where 


It  is  clear  that  the  eigenvalues  “i  and 
eigenfunctions  Ux  are  functionals  of  the  design 
variables  S,  and  they  depend  on  S  only  through 
the  damping  operator  C.  If  S  is  changed  by  an 
amount  4S,  then  u>x  and  Ux  will  change  by  a«i 
and  AUi,  respectively,  where 


■  wi(S  *  4S)  “  “i  (S) 
-  Ux  <S  +  55)  -  Ux  (S) 


(4) 


Assuming,  possesses  at  least  two  Gateaux 
derivatives  with  respect  to  S,  and  Ux  at  least 
one  Gateaux  derivative  with  respect  to  S,  then 
Eq.  (4)  may  be  expressed 

a«x  m  5«»x  ♦  4^x  ♦  0(4S3) 

2  <5 
4Ux  -  4Ux  ♦  0(5S  ) 

where  4*x  and  4Ux  are  linear  in  4S  and  is 

quadratic  in  4S.  A  major  objective  of  this 
paper  is  to  determine  explicit  representations 
for  5wt  ,  42u>x  and  4Ux  in  terms  of  the  eigen¬ 
values,  eigenfunctions  and  operators  appearing 
in  Eq.  (2a) . 


_  ,  .  ^  X,., 

,(G,  5CU.) 

(  Kx  ' 

-  4wx 

(vj 

\0  ' 

The  solution  to  Eq.  (9)  can  be  conveniently 
obtained  by  introducing  the  eigenvalue  equation 
which  is  adjoint  to  Eq.  (7) .  Toward  this  end, 
define  the  eigenvalue  problem 


where  (Xx  Yx )  is  the  adjoint  eigenvector.  It 
is  straightforward  to  show  that  Eqs  <7;  and 
(11)  admit  the  same  eigenvalues  ,  and  that 
the  adjoint  eigenfunctions  are  orthogonal  in 

the  sense  that 

•(:;  )  *  <*1'V«  *  -  0 

for  ii*j.  At  this  point,  we  are  at  liberty  to 
arbitrarily  select  two  normalization  factors, 
one  for  each  of  the  eigenvectors.  Consequent¬ 
ly,  without  loss  in  generality,  cr.e  set  of 
normalization  factors  will  be  selected  so  that 
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<Xi  ii) 


sij 


(12) 


where  is  the  Kronecker  delta. 

In  order  to  compute  the  first  variation 
4<*»i,  we  observe  that 


<*i  *i> 


0. 


(13) 


Equation  (13)  is  established  by  substituting 
Eq.  (9)  into  the  left  side  of  Eq.  (13),  and 
changing  the  order  of  integration  in  the 
resulting  double  integrals  using  the  symmetry 
properties  of  C,  M  and  G.  By  substituting  Eq 
(10)  into  Eq.  (13),  4“>i  is  immediately  deter¬ 


mined.  Thus 


«»i  “  “j2  <Xi,  <G,  4CU^ )  n  )fl  (14) 

Equation  (14)  provides  an  explicit  represen¬ 
tation  for  the  first  variation  of  the 
eigenvalue. 

The  first  variation  of  the  eigenfunctions 
are  not  so  simply  obtained.  By  virtue  of  Eqs. 
(10)  and  (14),  the  right  side  of  Eq.  (9)  is 
explicitly  known.  It  therefore  remains  to 
solve  the  inhomogeneous  linear  integral  aqua¬ 
tion  (9)  for  Vj^)7  .  Toward  this  end,  let 

the  solution  be  represented  in  terms  of  the 
eigenfunctions,  i.e. 


« 


(15) 


T 

Further,  let  (F^  Hi)  also  be  expanded  in  terms 
of  these  same  eigenfunctions.  Thus 


(16) 


Now  substituting  Eqs.  (15)  and  (16)  into  Eq. 
(9),  simplifying  the  result  using  Eq.  (6),  and 
comparing  term  by  term,  we  find  that 


U 


(17) 


The  coefficients  fij,  obtained  directly  from 
Eqs.  (12)  and  (16),  satisfy 

«ni  “  <*n  (18) 

Now,  substitution  of  Eqs.  (10)  and  (14)  into 
Eq.  (18)  immediately  yields 


U )  .  U) 

4U,  =  £  — — D  (x  ,  (G, 5CU  )  )  U 
1  r*i  V“n  n  1  n  a  n 


a  u 
i 1 


(20) 


Similarly,  it  can  be  shown  that  the  first 
variation  of  the  adjoint  eigenfunction  :X, 
satisfies 


U  U) 

4X  =  £  -i— -  ( X .  ,  (G,  6CUn)  )  X  -  a  X 

1  r*i  i~"n  1  n  n  n  n  n  i 


(21) 


The  constants  a^,  as  yet,  are  indeterminate. 
This  is  because  only  one  normalization  scheme 
has  been  involved  for  the  two  sets  of  eigen¬ 
functions.  Prescription  of  the  second  nor¬ 
malization  factor  will  render  aw  determinate. 

2  The  second  variation  of  the  eigenvalue 
4  can  be  computed  directly  from  Eqs.  (14), 
(20)  and  (21)  .  Although  some  manipulations  are 
required,  it  is  nevertheless  straightforward  to 
establish  that 


4\  =  2u.3(X-i,(G,4CUiV[l2  ♦  <V 


CUiVn 


+  2w'2  1  ^(Y<c'scWA,( 


n*i  i  n 


(22) 


It  is  interesting  to  note  that  neither  the 
first  nor  the  second  variation  in  the  eigen¬ 
value  depends  upon  the  particular  second  nor¬ 
malization  factors  selected  for  the  eigen¬ 
functions  . 

APPLICATION 

Equations  (14)  and  (22)  provide,  in  very 
general  form,  explicit  expressions  for  the 
first  and  second  variations  of  the  complex 
frequency  for  unforced  viscously-damped 
vibrations  of  elastic  systems.  Specific  usage 
of  the  equations  will  vary  with  the  applica¬ 
tion.  In  this  section,  we  shall  be  concerned 
with  the  development  of  a  closed  form  approx¬ 
imation  for  for  structures  in  the  presence 
of  "moderate"  damping. 

Without  loss  in  generality,  we  may 
restrict  our  treatment  of  damped  elastic 
systems  to  consideration  of  the  quadratic 
eigenvalue  equation 


f 


ni 


’“i 

o 


n¥i 

n-i 


(19) 


The  first  variation  40* ,  obtained  from  Eqs. 
(15),  (17)  and  (19),  becomes 


Ei  -  G  C  oi  -  uii2  G  2i  (23) 

where  G  is  a  diagonal  matrix  specified  in  terms 
of  the  system' s  natural  frequencies  .... .  Speci¬ 
fically, 


ij 


-2 

*i  41J 


(no  sum) 


(24) 
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Since  there  are  an  infinite  number  of  natural 
frequencies  associated  with  vibrating  distri¬ 
buted  parameter  structures,  G  will  also  contain 
an  infinite  number  of  elements.  Likewise,  each 
eigenvector  is  an  infinitely  dimensional 
column  vector.  The  matrix  Q  has  elements 

Cij  -  (Wi,  C  Wj)n  (25) 


Clearly  C  (1)  -  C.  Let  ^(S)  and  ;x(S)  be, 
respectively,  the  eigenvectors  and  eigenvalues 
associated  with  £<S) .  Now 

S^S)  -  ^(0)  +  Si'(0)S  +  +  ... 

where  C^l)  is  the  desired  eigenvalue  *ix. 
Therefore 


where  Wx  is  the  normalized  mode  shape  corre¬ 
sponding  to  the  (undamped)  natural  frequency 

°i* 

In  order  to  establish  the  reduction  of  Eq. 
(3)  to  the  form  (23) ,  it  is  necessary  to  recall 
the  eigenvalue  equation 


Wn  * 


(G, 


M  V0 


(26) 


and  the  orthonormal  relations 


(WP,  M  Wn)n  -  < ^  (27) 

for  undamped  vibrating  elastic  structures.  It 
is  assumed  that  nn  and  wn  are  known  for  the 
structure  under  consideration.  Since  the  set 
of  eigenfunctions  (Hn)  form  a  complete  basis, 
the  complex  eigenfunction  may  be  expressed 


®i 


r  a. 
n  in 


(28) 


where  each  is  a  complex-valued  constant. 

Now  substitute  Eq.  (28)  into  Eq.  (3),  take  the 
I>2  inn**-  product  of  each  side  of  the  resulting 
equation  with  M  wp,  change  the  order  of  the 
double  integrals  that  arise  and  simplify  the 
result  using  Eqs.  (26)  and  (27).  The  resulting 
equation  is 

“ip"  '  “i  £  £°in  <V  C  Vo*  “i  "pNp  <29> 


“i  “  “i(°)  +  “i’(0)  +  T  Sx"(0)  + 


(32) 


To  compute  <*^(0) ,  set  C  -  C0  in  Eq.  (23), 

thereby  uncoupling  the  equations.  The  solu¬ 
tions,  obtained  immediately,  are 


“i<0)  - 


'ii 


-40 


(33) 


Alternatively,  S^(0)  are  the  eigenvalues 
obtained  by  ignoring  the  off-diagonal  elements 
in  the  damping  matrix  C.  Indeed,  Hutton  (1981) 
suggests  the  approximation 


“x  *  fit(0)  (34) 

may  be  used  whenever  the  off-diagonal  elements 
of  C  are  sufficiently  small.  Finally,  for 
future  reference,  it  is  noted  that  the  eigen¬ 
vector  _8X(0)  associated  with  £^0)  has 
components  which  may  be  taken  as  the  Euclidean 
unit  base  vector  ,  so  that 

“ip(0>  -  «ip  ,  ot( 0 )  -  ex  (35) 

It  now  remains  to  calculate  the  sensi¬ 
tivity  derivatives  appearing  in  Eq. (32) .  The 
first  derivative  u’(0)  is  computed  directly 
from  Eq.  (14)  by  observing 

«“i-  £^(0)  AS  .  5C  ■  4t  SS  (36) 


which  is  precisely  Eq.  (23)  expressed  in  compo¬ 
nent  form. 

It  should  be  pointed  out  that  with  similar 
arguments  Eq.  (23)  can  also  be  established  for 
discrete  finite-dimensional  eigenvalue  problems 
as  well.  In  this  case  the  l>2-inner  product  is 
replaced  everywhere  by  the  usual  dot  product  of 
complex  vectors. 

Ne  now  take  up  the  development  of  an 
approximate  solution  to  Eq.  (23)  using  the 
eigensensitivity  theory  developed  in  the 
previous  section  of  this  paper.  Toward  this 
end,  it  is  convenient  to  define  a  one  parameter 
family  of  damping  matrices  C  (S)  where 

C  (S)  -  CD  ♦  S  AC  (30) 


in  which  CD  contains  only  the  diagonal  elements 
of  C  and  a?  contains  only  the  off-diagonal  ele¬ 
ments  of  C.  Thus 

CD  -  C  xj  s  ij  (no  sum) 


and 


ACXj 


Cij 

0  i=j 


(31) 


Further,  a  comparison  of  Eqs.  (23)  with  Eqs. 
(3),  (6)  and  (11)  reveals  that  M-l  and 


Oi  "  ii  •  “  “i  G  C  £x  ~  “x  9  ii  (37a) 

Vi  -  “x  ii  (37b) 

xi  *  ’  “i  ?  ?  1  G  ^  (37c) 

*i  -  -  “i  <3  ii  ( 37d) 

Since  C(0)„  •  5°  i*  a  diagonal  matrix, 

Q.  C(0)  -  C (0) Q,  and  Eqns.  (35)  and  (37a, c) 
imply  that 

xi  (0)  “  ui  ®i  <30) 

where  ux  is  a  complex-valued  scalar.  It  can 
also  be  established  directly  from  Eqs.  (35)  and 
(37b, d)  that 


(0)  -  -  «x  sx  <°>  ? 

Vx  (0)  -  -  “x<°)  «1 


(39) 
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The  scalar  Mi  can  be  computed  from  the  ortho¬ 
normality  condition  (12)  and  Eqs.  (35),  (38), 
and  (39) .  Thus 

«i  -  0i2/(sJi2  -  ^(0) )  (40) 

After  substituting  Eqs.  (24),  (36)  and 
(37a, c)  into  Eq.  (14) 

V<°>  -  5?(0)  *1  •  ?  «i  (41) 

results.  After  carrying  out  the  indicated 
operations  in  Eq.  (41),  we  obtain 

«i<0)  -  ni  S2(0)(j-.2  4  ,  Q<  (42) 

Equation  (42)  validates  Hutton's  (1981) 
approximation  (34)  whenever  <j£  is  very  small, 
for,  in  this  case,  the  error  in  Eq.  (34)  is  of 
order  ( AC) 2. 

It  still  remains  to  calculate  u"(0)  from 
Eq.  (22) .  He  observe  that 

«i”  (0)  «s,2 


n-i,  it  must  be  recognized  that  there  are  two 
values  for<%  (0)  corresponding  to  each  n  which 
must  be  included  in  the  summation. 

Also,  it  should  be  noted  that  the 
sensitivity  derivatives  are  evaluated  at  S«0 . 
These  derivatives  exist  only  for  distinct 
zeroth  order  roots  u>  £0) ,  whether  complex  or 
real-valued.  And  finally,  it  can  be  easily 
shown  that  the  other  terms  appearing  in  the 
denominator  of  the  expression  for  (0), 
namely  £  j(0)  never  vanish. 

EXAMPLES 

(a)  As  a  very  simple  illustration,  consider 
the  elementary  non-dimensional  quadratic 
eigenvalue  equation  corresponding  to  the  two- 
degree  of  freedom  system  in  which  C  «  l, 

C  22  -  2,  C  12  -  -0.75,  (1?  -  2,  -  4.  The 

zeroth  order  approximation,  obtained  by 


ignoring  C  12,  is 

£  ^O)  -  -0.5  t 

1.9365  i 

(i2(0)  -  -1.0  i 

1.0000  i 

and 

-  0 


On  the  other  hand,  the  frequencies  calculated 
from  the  second  order  approximation  (45)  are 


Now,  we  combine  Eqs.  (22),  (24),  (33),  (35), 
(36),  (37a, c),  (40),  (41)  and  (43)  to  obtain 


,  ,  fiT(O)  r  (AC  r 

i  a»(oj  =  1  .  r  - dl_ 

°i  '“V0*  n?i  L  <^(0)-^ 


“1  «  -0.3313  1  1.8639  i 

(47) 

“2  *  -1.1688  *  0.9438  i 

By  way  of  comparison,  the  exact  eigenvalues  are 
“1  -  -0.3434  i  1.8863  i 


(44) 


Finally,  after  substituting  Eqs.  (42)  and  (44) 
into  Eq.  (32),  we  obtain  the  result 


»l(0)  ♦ 


<*C_,)' 


°f  ->)  nrfi  fi^W-WnfO)  *  ~ 


S  S(0> 


(45) 

where  <2^0)  is  given  by  Eq.  (33)  .  If  desired, 
the  approximation  (45)  can  be  made  more 
accurate  by  including  still  higher  order 
derivatives  of  &i(S).  In  this  case,  J3  can 
be  calculated  directly  from  Eqs.  (20-22). 
Although  the  procedure  is  straightforward,  the 
final  expression  for  i^»iis  cumbersome  and 
consequently  will  be  omitted  here. 

A  few  observations  regarding  the  usage  of 
Eq.  (45)  is  in  order.  The  frequencies  «i(0) 
always  occur  in  pairs:  if  they  are  complex¬ 
valued,  then  the  pairs  are  complex  conjugates; 
but  if  they  are  real  valued,  then  the  pairs 
generally  consist  of  distinct  negative  numbers. 
When  summing  over  all  values  of  n,  other  than 


"2  -  -1.1576  t  0.9146  i 

Obviously,  the  second  order  approximation  (47) 
is  far  superior  to  the  zeroth  order 
approximation  (46) .  As  a  general  rule,  the 
most  important  value  is  the  smallest  (in 
magnitude)  real  component  of  u>if  in  this  case 
Re  (u>i).  The  zeroth  order  approximation  for 
Re  (u^)  is  in  error  by  50.9%,  while  the  second 
order  approximation  for  Re  (w.)  is  in  error  by 
3.6%. 

The  accuracy  of  Eq.  (47)  is  somewhat 
remarkable  considering  that  C12  is  not  small 
compared  with  either  of  the  diagonal  elements. 
However,  it  can  be  shown  that  for  systems  with 
two  degrees  of  freedom,  Eq.  (45)  is  accurate  to 
within  terms  of  order  (AC)4,  and  hence  the 
rapid  convergence  of  the'series  approximation 
(32)  . 

(b)  A  somewhat  more  involved  quadratic 
eigenvalue  problem  is  provided  by  the  four 
degree  of  freedom  system  in  which  C 1  .  «  2, 

C  22  ■  4,  C33  •  1,  C44  *3,  Ctj  -  -0.25  for 
i  t  j,  end  0  J -  oj  -  2,  -  n  3 .  In  this 

example,  the  off  diagonal  damping  elements  are 
small,  but  not  negligible,  compared  to  the 
diagonal  damping  coefficients.  There  are  now 
four  pairs  of  solutions  ui  j,  of  which  one  pair 
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ia  real-valued  and  tha  others  ara  complex- 
valuad.  Tha  zeroth  order  solution,  obtained 
immediately  from  Eq.  (33),  is 


U»  J  (0) 

m  2  (0) 

-  -1  i  i 

-  -0.5858, 

-3.4142 

i 

(48) 

w3  (0) 

-  -0.5  * 

1.6583 

i4(0) 

-  -1.5  * 

0.8660 

i 

However,  using  the  second  order  approximation, 
tha  solution  is 


"l 

ss 

-1.0094 

t 

1 . 0656i 

"2 

» 

-0.5718 

9 

-3.5332 

“3 

a 

-0.4493 

1 

1 . 6605i 

"4 

« 

-1.4888 

1 

0.76161 

It  can  be  shown,  although  with  considerably 
more  effort  than  was  required  in  tha  two-degree 
of  freedom  example,  that  the  exact  frequencies 
are 


-1.0142 

t 

1.0506 

i 

*2 

-0.5744 

9 

-3.5062 

(50) 

U»3  m 

-0.4377 

1 

1.6584 

i 

U <4  - 

-1.5077 

t 

0.7591 

i 

Once  again,  the  largest  error  in  the  zeroth 
order  frequency  occurs  in  the  least  damped 
frequency  -  in  this  case.  Re  (W3}.  Here,  the 
quadratic  approximation  produces  an  error  of 
2.7%  compared  with  14.2%  for  the  zeroth  order 
approximation. 

(c)  As  an  illustration  of  a  continuous  system, 
consider  the  unforced  vibrations  of  a  uniform 
cantilever  beam  subject  to  distributed  damping. 
Specifically,  the  quadratic  eigenvalue  equation 
is 

EIOIV(x)  +  C(x)  <u  0(x)  +  P^Of*)  -  0  (51a) 

with  boundary  conditions 

U  (0)  -  U  '(0)  -  0  "(L)  -  U  (L)  -  0  (51b) 

It  will  be  assumed  that  uniform  damping  occurs 
over  the  outer  half  of  the  beam,  so  that 

C  (x)  -  C  0  H  (x-L/2)  (52) 

where  C0  is  the  damping  coefficient  and  H  is 
the  Heaviside  function.  Further,  El  is  the 
beam's  uniform  flexural  stiffness,  L  is  its 
length,  p  is  its  mass  density  per  unit  length 
and  U  is  the  response  eigenfunction.  It  is 
convenient  to  cast  Eq.  <51a)  into  non-dimen¬ 
sional  form  by  setting 

Z  •  U/L  ,  y-  x/L,  a2  •  1, 4  (pu,VEI) 

S  •  C/2  Vo"EI  ,  S„-  C,/2  vSTl  (53) 

so  that  Eqs .  (51a, b)  and  Eq.  (52)  become 


z Iv  +  2  paz  +  a2  z  -  0 

Z(0)  -  Z *  (0)  -  Z"(l)  -  Z'M  (1)  -  0  (54) 

0  -  0.  H  (y  -  1/2) 

Here,  a  is  the  desired  eigenvalue,  Z  the 
corresponding  eigenfunction,  and  2s  is  the 
damping  operator. 

In  order  to  use  the  quadratic  approx¬ 
imation  (45),  Eq.  (54)  must  be  cast  into  the 
form  (23) .  It  is  therefore  necessary  to  first 
obtain  the  natural  frequencies  fln  and  corre¬ 
sponding  mode  shapes  Wn(y)  for  the  undamped 
system,  i.e.  0  ■  0.  These  results,  however, 
are  readily  available.  If  we  let 

,,ij  *  Vy)dy  (55) 

the  damping  matrix  becomes 

cij  *  2  B„  ntj  (56) 

Direct  computations  of  Sj,  provide  specific 
values  of  C^j.  The  zerotn  order  frequencies 
now  become 

OjlO)  -  *  J,  «ltii  i^n  2  -  MnJ  (57) 

while  the  second  order  approximation  to  the 
frequencies  can  be  obtained  from  Eq.  (44)  using 
specific  values  for  e»  . 

Results  for  the  first  six  eigenvalues  are 
presented  in  Table  1  for  light  damping  B„  -  1.5 
and  moderate  damping  P.  -  3.0.  For  B  „  «  1.5, 
the  zeroth  order  frequencies  (57)  are  very 
close  to  the  more  accurate  second  order 
approximation,  denoted  by  a}  .  The  maximum 
difference  in  the  real  parts  of  the  frequencies 
occurs  in  CI3,  the  least  damped  mode,  and  is 
approximately  4%. 

F<jr  9*  -  3.,  the  difference  between  a  t(0) 
and  a  t  is  more  pronounced.  As  shown  in  Table 
2,  both  eigenvalues  corresponding  to  are 
real,  while  the  remaining  ones  are  complex¬ 
valued.  Except  for  the  second  and  third 
frequencies,  the  quadratic  approximation  a  * 
adds  relatively  little  to  the  much  simpler 
solution  ai(0).  But  the  Re  {  a  3(0 )  )  predicts 
16%  more  damping  than  Re  (  a3* ) ,  while  the 
Re  ( 3?(0)  )  predicts  10%  less  damping  than 
Re  (a;  ). 

in  the  absence  of  an  exact  solution  to  Eq. 
(54),  it  is  difficult  to  precisely  determine 
the  error  introduced  by  truncating  the  series 
(32)  after  quadratic  terms.  The  exact  solution 
is  the  limiting  case,  as  N  —  *,  of  an  N  degree 
of  freedom  system  obtained  by  considering  the 
sub-matrices  of  C  and  G  in  which  i,  N. 

Unfortunately,  it  is  difficult  to  solve, 
exactly,  the  quadratic  eigenvalue  problem  when 
N  is  not  a  small  number.  Nevertheless,  by 
examining  the  exact  eigenvalues  for  the  finite 
dimensional  problem  for  small  values  of  N,  an 
assessment  of  the  relative  accuracy  of  a,*  can 
be  made. 
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fli 

<0) 

a 

* 

i 

1 

-1.424 

r  1.2201 

-1.441 

1 

1.2161 

2 

-0.838 

$  4.6191 

-0.858 

♦ 

4 . 684i 

3 

-0.739 

t  7 . 820i 

-0.712 

♦ 

7 . 716i 

4 

-0.748 

*  10.9701 

-0.752 

♦ 

10.9961 

5 

-0.750 

t  14 . 117i 

-0.744 

t 

14 . 058i 

6 

-0.750 

t  17 . 262i 

-0.752 

♦ 

17 . 279i 

7_ 

-iM.ion 

_ r<?,7j7 

♦ 
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S' 

i 

a 

1  1 

10) 

a 

1 

1 

-4.991 

» 

-0.704 

-5.187, 

-0.699 

2 

-1.676 

t 

4 . 385i 

-1.867 

♦ 

4 . 623i 

3 

-1.479 

t 

7 . 714i 

-1.273 

t 

7.3411 

4 

-1.497 

1 

10.8931 

-1.526 

t 

10.9941 

5 

-1.501 

1 

14.057i 

-1.449 

t 

13 . 824i 

6 

-1.500 

1 

17 .2131 

-1.511 

t 

17 . 277i 

7_ 

-C1..-S.W 

♦ 

20-3651 

_ sLIZS 

t 

iQ-lQU 

N- 

2 

N-3 

N-< 

t 

N-5 

N-6 

n  -5 

.137, 

-5. 

.147 

-5. 

150 

-5. 

.  151 

-5. 

152 

a  1  -0 

.700 

-0. 

.699 

-0. 

699 

-0. 

699 

-0. 

699 

Re(a2)  - 

1.605 

-1. 

.934 

-1, 

.925 

-1. 

.942 

-1 . 

942 

Im(a2» 

4.355 

4. 

.659 

4  . 

,661 

4. 

,673 

4  . 

,673 

Re  { a  3) 

_ 

-1 

.145 

-1, 

.208 

-1. 

.  188 

-1  , 

.  190 

Im<03) 

- 

7, 

.132 

7  , 

.273 

7  . 

.258 

7  , 

.261 

Re  { a  4) 

- 

- 

-1 , 

.442 

-1 

.569 

-1  . 

.559 

Xm  ( a  4} 

- 

- 

10 

.678 

11 

.047 

11 

.033 

Re(a,l 

• 

_ 

-1 

.376 

-1 

.429 

im[a\) 

- 

- 

- 

13 

.566 

13 

.782 

Re{a  6> 

- 

- 

- 

- 

-1 

.455 

Ie(a  ft) 

- 

- 

- 

- 

16 
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8UMKART 

We 

have  presented 

an 

approximate 

solution 

to  the  quadratic  eigenvalue  problem  by  expand¬ 
ing  the  eigenvalues  in  a  series  about  A  c-0  and 
truncating  after  quadratic  terms.  It  was 
shown,  by  examples,  that  whenever  AC  is  small 
but  not  negligible,  the  series  solution  (45) 
provides  reasonable  estimates  of  the  true 
eigenvalues  even  when  the  solution  obtai  »d  by 
ignoring  AC  is  substantially  in  error. 
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ABSTRACT 

A  new  closed-form  approximation  for  the 
fundamental  frequency  of  symmetric 
rectangular  laminates  subject  to  all 
combinations  of  hinged  and  clamped 
boundary  conditions  is  presented.  The 
distributed  parameter  eigenvalue  equation 
is  cast  in  an  equivalent  infinitely- 
dimensional  discrete  form.  The  stiffness 
and  mass  matrices  are  each  decomposed  into 
the  sum  of  two  matrices,  one  of  which  is 
diagonal  while  the  other  contains  zero 
diagonal  elements.  Design  sensitivity 
analysis  is  used  to  expand  the  desired 
eigen frequency  in  a  Maclaurin  series  of 
the  zero  diagonal  matrices.  The  general 
formula  thus  obtained  is  then  specialised 
to  rectangular  symmetric  laminated  plates. 
The  remarkable  occuracy  of  this  new 
formula  is  established  by  numerical 
comparisons  of  results  calculated  from  it 
to  those  obtained  from  the  conventional 
Rayleigh-Rits  method. 

INTRODUCTION 

Methods  to  accurately  determine  the 
frequencies  of  laminated  plates  is 
necessary  for  both  analysis  and  design  of 
symmetric  laminated  composite  plates. 

While  exact  closed  fora  solutions  for  the 
natural  frequencies  are  available  only  for 
simply-supported  orthotropic  laminates 
(1],  very  accurate  approximate  solutions 
are  also  available  for  orthotropic 
laminates  meeting  other  boundary 
conditions  (2]. 

For  more  general  anisotropic  plates,  the 
literature  contains  several  closed  form 
approximations  for  the  fundamental 
frequency.  Bert  obtained  simple  formula 
for  symmetric  rectangular  laminates  for 


simply  supported  (3)  and  clampea  edges  [ 4 } 
by  assuming  a  linear  relationship  between 
the  eigenvalue  and  the  flexural 
stiffnesses.  A  different  approach  was 
adopted  by  Reiss  and  coworkers.  Reiss, 
Ramachandran  and  Qian  (5J.  using  exact 
sensitivity  derivatives,  obtained  an 
approximate  solution  for  the  natural 
frequencies  of  a  four-layer  symmetric 
simply- supported  angle-ply  laminate.  The 
frequencies  obtained  were  presented  in  the 
form  of  a  truncated  Fourier  series  in  the 
ply-angle.  Reiss  and  Qian  [6,7]  then 
generalized  those  results  to  include  all 
symmetric  laminates,  none  of  whose  edges 
were  free. 

In  the  present  study,  a  new  and  much  more 
general  eigensensitivity-based  method  for 
obtaining  approximate  formula  for  the 
frequency.. response  of  structures  is 
employed.  The  method,  recently  developed 
by  Reiss,  Qian  and  Aung  (8)  for  damped 
systems  and  Reiss  and  Qian  [9]  for  self- 
adjoint  systems  requires  recasting  the 
distributed  parameter  equation  into  an 
equivalent  discrete  set  of  algebraic 
equations.  The  new  formula  for  the 
fundamental  frequency  of  symmetric 
laminates  is  obtained  by  specializing  the 
general  result  [9]  to  the  specific  problem 
under  consideration. 

REVIEW  OF  FREVIOOS  R1SULTS 

At  this  point,  it  is  useful  to  review  the 
general  methodology  previously  developed 
by  Reiss  and  Qian  (9). 

The  general  eigenvalue  problem  under 
consideration  is  representable  in  the  form 
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T'£ Tv,  «  l,^  in  Q  ( 1  > 


subject  to  appropriate  Mixed  boundary 
condition*.  Rare  T  and  T*  are  L,  (0) 
adjoint  differential  operator*,  E  and  M 
are,  respectively,  the  stiffness  and  mass 
operators  which  are  each  positive  and 
symmetric  with  respect  to  the  L,  (O)  inner 
product.  Also  X,  is  the  eigenvalue 
associated  with  the  eigenfunction  u(. 

While  there  are  many  ways  in  which 
Equation  (1)  may  be  converted  into  a 
discrete  fora,  the  Raleigh-Ritz  method  is 
particular  useful  in  view  of  the 
objectives  of  this  paper.  Thus,  let  (V») 
denote  a  complete  set  of  kinematically 
admissible  basis  functions,  so  that 


ui  •  E  <2> 

fl*l 


Stationarily  of  the  Rayleigh  quotient 
requires 


where  the  square,  infinitely  dimensional 
autrices  (K)  and  (M]  contain  elements 


°,f  =<lu*tion  (5).  it  Is  evident  the 
VI  1,1  ar®  diagonal  matrices,  while 

[kj  and  fm)  contain  only  zero  values  on 
the  diagonal.  It  is  convenient  to 
consider  Xt  to  be  a  scalar-valued  function 
of  the  arguments  I k J  and  [«J .  Indeed,  if 
is  expanded  as  a  Haclaurin  series  in 
these  arguments,  it  can  be  shown  that  X. 
satisfies  19)  1 


*i  -  n;}  -  ui!  y  *m  -  Ha 

nJr  *-»  - 

(7) 


ord^T  r#t?ln*d  series .  Higher 

th*  ®xP*n3ion  may  be  9 
eadily  calculated;  however,  only  bv 
compromising  the  simplicity -of  Equation 

SYMMETRIC  LAMINATED  PLATES 

The  eigenvalue  equation  for  a  freelv 
vibrating  symmetric  laminate  is  Y 


°»iu«.  ♦  O,, <W  »  2  (O,,  »  2 u^, 
4  i>,fUJoar  *4  Di%  uw  *  p  u*  u 


(R) 


KtJ  • 

«  (fj.Mfj) 


(*) 


In  Equation  (4),  (  •,  )  denotes  the  L,(fl) 
inner  product.  Further  a*.  is  the  n-th 
component  of  the  infinitely  dimensional 
column  vector  (a,). 

By  introducing  new  matrices  (K*) ,  (kj,  [M*J 
and  (*)  defined  by 


where  u  denotes  the  mode  shape,  o>  is  the 
corresponding  frequency,  p  is  the  mass 
density,  Dl#  are  the  flexural  stiffnesses 
11],  x  and  y  are  the  coordinates  in  the 
plane  of  the  laminate,  and  indices 
following  u  denote  differentiation  with 
respect  to  the  indicated  argument. 

In  order  to  apply  Equation  (7),  it  is 
necessary  to  cast  Equation  (8)  in  the  form 
of  Equation  (1),  that  is,  the  operators  T, 
T  ,  E,  M  and  the  basis  functions  (y,)  must 
be  identified.  It  is  straight-forward  to 
show  that 


K%*  -  K„  84,  (no  sum) 

Mjj*  -  Mtj  St,  (no  sum) 

k»l  ■  K»»  -  K»»  •»>  <“«  *“> 


(5) 


T 


a*/a** 

avay® 

2d*/dxdy 


T*  -  Tf,  M  -  p. 


X  -  (J  (9) 


*U  -  M„  -  M,,  •„  (no  sum) 

where  8„  is  the  Kronecker  delta,  it 
follows  that  Equation  (3)  becomes 

((A®)  *  IJt) ) la,)  •  IjUM*)  ♦ 

(6) 


and,  since  no  edges  are  free,  the  stiff¬ 
ness  operator  may  be  selected  as 


On 

0. 

Du 

E 

m 

0. 

0„ 

D„ 

o,» 

D„ 

0 
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where 


where 


D,  •  0U  ♦  2  0M 

C.1  - 

(X.’.Xj) 

C,,  -  (Y,  ,  Y,') 

As  is  common  in  plat*  theory, 
functions  will  be  selected  as 

the  basis 
products  of 

G,,  - 

9i»  -  (7i‘.Y,) 

the  beam  shape  functions,  !•* 

• 

H„  - 

(X,  ".X,  ) 

h„  -  (Yj  ",  Yj‘) 

♦-»  «*.y>  •  9* .  <*/•>  (ll) 


whtr«  th*  plat*  domain  Is  0  S  *  £  •# 

0  S  y  S  b  and  R  -  a/b  i»  the  laminate'  a 
aspect  ratio.  Th*  baa*  shape  function a 
satisfy 


x r  -  *\  *.  *  0  0  a  x/a  a  1 
y/*  -  »*  r,  •  0  0  a  y/b  a  1 


In  view  of  Equations  (15-16),  tha  terms 
appearing  in  Equation  (14)  are 

pa’K,,,,  -  Du|i,#  ♦  R*OtlU,4  ♦  2D,R»CUCU 

^llll  *  Kn,  «  1 

P*'*"5—,  -  ♦  R’DjjU,*  ♦  2DtR*C-,CMl 

Pa‘*«H  "  2D^*CBlc,,  «•  2RD„glll  (H., -!»,„) 

♦  2R,0„G.l  (hu-hBl) 

<Vn  -  0  (17) 


uhara  (  )"  denotes  the  fourth  derivative 

with  respect  to  the  indicated  non-dimen¬ 
sional  argument.  It  is  worth  noting  that 
v  t  y  $  ^  and  k-nown  for  tftch  of 

boundary  conditions.  rurther,  if  th* 
basis  functions  are  normalised  with 
respect  to  th*  mass  density  p,  it  follows 
that  th*  shape  functions  eust  satisfy 


rinally,  substitution  of  Equations  (17) 
into  (14)  produces  the  desired  result 

pa‘XM  -  pa'to,,*  -  D, i|St*  ♦  R‘0,,0,*  ♦  aD^C.jC,, 

'  4R\  E  L°*c»c»i  *  (H.t-Hu) 

t) 

♦  0,^*6.,  (h„-h„)  ]  */  toM  (p,*-p,‘) 


(Xt.X,)  ■  OVI-jl  «  P  l/i  «„  <12> 


♦  D„R‘ («„•-«,•)  ♦  ZO.R'fC^.-CuC,,)! 


It  remains  to  determine  th*  mass  and 
stiffness  elements  from  Equation  (4) . 

Since  th*  basis  functions  in  Eq.  (11)  are 

double  subscripted,  it  follows  that  U* 

proper  form  for  Eq.  (4)  is  DXSC039I0B 


3  (13) 

while  Equation  (7)  becomes 

.  ^Iim  .  1  re 

11  *1111  Mfm 

(14) 

Substitution  of  Equations  (8-U)  into  (13) 
readily  yields 


Equation  (IS)  provides  the  desired 
approximation  for  the  fundamental 
frequency  of  laminated  rectangular  plates 
subject  to  all  combinations  and  clamped 
and  simply-supported  boundary  conditions. 
In  the  derivations  of  this  equation,  it 
was  implicitly  assumed  that  the  natural 
frequency  cs^t  is  in  fact  the  fundamental 
frequency.  It  can  be  established  that  for 
the  boundary  conditions  under  consider¬ 
ation,  this  assumption  is  indeed  correct, 
rinally,  it  should  be  pointed  out  that  th* 
matrices  (C),  (e),  (0).  (g),  (I)  and  (b) 
and  beam  frequencies  p*  and  ut  are  boundary 
condition  dependent;  hence  these  values 
account  for  the  effect  of  th*  boundary 
conditions  in  Equation  (18) . 


*!/■  ■  ®i»  fl*) 

p  a4  X,jm  ■  (0|iR*  ♦  ®40u  7j) 

♦  2^  **  Cyp*  ♦  2 

*  21'Dh  Gm  (16) 


Apart  from  th*  term  proportional  to  D,  in 
th*  numerator  of  th*  double  sum  term 
appearing  in  Equation  (18),  th*  same 
expression  for  cs^  was  obtained  by  Reiss 
and  Qian  (7] .  In  their  approach,  tt,t  was 
treated  as  a  function  of  D,»  and  D>(,  and 
th*  sensitivity  derivatives  of  (0^  with 
respect  to  these  arguments  were  used  to 
generate  th*  approximate  solution. 
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However,  sines  thsss  derivatives  could  not 
be  exactly  evaluated,  a  somewhat  hueristic 
approach  was  used  to  obtain  the  ultimate 
result.  In  the  approach  of  this  paper, 
however,  exact  sensitivity  derivatives 
have  been  obtained,  and  consequently 
Equation  (18)  is  exact  to  within  the  order 
of  terms  retained  in  the  Maclaurin  series. 


finally  the  edge  y"b.  The  columns 
labelled  k,  are  calculated  using  the  tut* 
method  with  N»9,  while  the  column  labelled 
kt  is  determined  from  Equation  (18)  also 
with  N«9.  Finally,  since  the  frequency  is 
symmetric  in  0,  results  are  limited  to  ply 
angles  between  0*  and  90*. 
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aspect  ratios  and  boundary  conditions. 
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TABLE  I  RESULTS  FOR  THE  SSSS  PLATE 


R 

*1 

kl 

*1 

kl 

B-E 

C-E 

0 

l 

38.895 

38.895 

37.954 

37.954 

is 

1 

41.010 

41.075 

39.815 

39.991 

30 

l 

45.509 

45.782 

44.338 

44.957 

4$ 

1 

47.789 

48.136 

46.784 

47.524 

60 

l 

45.509 

45.  782 

44.338 

44.957 

73 

1 

41.010 

41.075 

39.815 

39.991 

90 

1 

38.895 

38.895 

37.954 

37.954 

0 

2 

15.590 

15.590 

11.687 

11.687 

15 

2 

18.117 

18.129 

14.945 

15.023 

36 

2 

23.139 

23.215 

21.236 

21.492 

65 

2 

27.784 

27.945 

27.023 

27.384 

60 

2 

31.560 

31.688 

31.803 

32.042 

75 

2 

34.362 

34.383 

35.341 

35.382 

90 

2 

35.451 

35.451 

36.707 

36.707 

0 

5 

11.389 

It. 389 

6.208 

6.208 

15 

5 

12.151 

12.151 

7.859 

7.866 

30 

5 

15.178 

15.187 

12.719 

12.778 

65 

5 

20.831 

20.859 

20.242 

20.307 

60 

5 

27.540 

27.557 

28.244 

28.272 

75 

5 

32.803 

32.806 

34.234 

34.339 

90 

5 

34.776 

34  .  776 

36.442 

36.442 

TABLE  II  RESULTS  FOR  1>C  CCCC  PLATE 

B-E 

C-E 

1° 

R 

*1 

*1 

*1 

0 

1 

84.172 

84.172 

84.284 

84.284 

15 

l 

82.260 

82.181 

81.592 

81.467 

30 

1 

78.658 

78.355 

76.642 

76.416 

45 

l 

77.014 

76.580 

74.474 

74.216 

60 

1 

78.658 

78.355 

76.642 

76.416 

75 

1 

82.260 

82.181 

81.592 

81.467 

90 

1 

84.172 

84.172 

84.284 

84.284 

0 

2 

3T  817 

32.816 

25.050 

25.050 

15 

2 

34.273 

34.241 

26.997 

26.942 

30 

2 

39.773 

39.621 

34.846 

34.704 

45 

2 

50.513 

50.328 

48.936 

48.782 

60 

2 

64.002 

63.861 

65.210 

65.050 

75 

2 

75.055 

75.010 

77.963 

77.898 

90 

2 

79.275 

79.275 

82.756 

82.756 

0 

5 

25.276 

25.276 

13.644 

13.644 

15 

5 

25.892 

25.894 

15.603 

15.590 

30 

5 

31.498 

31.478 

25.746 

25.688 

45 

5 

44.918 

44.872 

43.749 

43.682 

60 

5 

61.212 

61.178 

62.955 

62.912 

75 

5 

73.911 

73.901 

77.268 

77.256 

90 

5 

78.651 

78.652 

82.536 

82.536 
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TABLE  111  RESULTS  FDR  THE  CSCS  PLATE 


B-E 

C-E 

R 

*1 

kl 

*1 

o 

l 

80.824 

80.824 

83.327 

83.327 

15 

l 

78.613 

78.637 

30.628 

80.686 

30 

1 

74.281 

74.276 

75.494 

75.500 

4.5 

l 

67.609 

67.638 

67.565 

67.621 

60 

1 

57.428 

57.481 

55.133 

55.236 

75 

1 

48.131 

48.134 

43.346 

43.318 

90 

1 

45.333 

45.333 

39.957 

39.957 

0 

2 

23.778 

23.778 

22.007 

22.007 

15 

2 

25.280 

25.275 

23.835 

23.833 

30 

2 

29.055 

29.077 

28.423 

28.461 

4.5 

2 

31.804 

31.834 

31.606 

31.653 

60 

2 

33.402 

33.348 

33.547 

33.477 

75 

2 

35.097 

35.051 

35.704 

35.625 

90 

2 

36.027 

36.027 

36.879 

36.880 

0 

5 

11.806 

11.806 

6.935 

6.935 

15 

5 

12.636 

12.641 

8.709 

8.716 

30 

5 

15.679 

15.682 

13.445 

13.435 

45 

5 

21.121 

21.077 

20.573 

20.505 

60 

5 

27.647 

27.599 

28.344 

28.279 

75 

5 

32.842 

32.826 

34.255 

34.234 

90 

5 

34.803 

34.804 

36.451 

36.452 

TABLE  IV 

RESULTS  FOR  TVE  SCCC  PLATE 

B-E 

C-E 

R 

*1 

S 

*1 

S 

0 

1 

61.839 

61.839 

59.340 

59.340 

15 

1 

62.322 

62.255 

59.334 

59.221 

30 

1 

64.755 

64.535 

61.623 

61.515 

45 

1 

69.466 

69.254 

67.296 

67.322 

60 

1 

75.034 

74.965 

74.264 

74.380 

75 

1 

80.008 

79.985 

80.820 

80.813 

90 

1 

82.104 

82.105 

83.701 

83.702 

0 

2 

29.487 

29.487 

20.052 

20.052 

15 

2 

31.377 

31.353 

22.997 

22.964 

30 

2 

37.934 

37.864 

32.769 

32.757 

45 

2 

49.665 

49.645 

48.174 

48.242 

60 

2 

63.6S8 

63.652 

64.99t 

65.011 

75 

2 

74.871 

74.859 

77.895 

77.883 

90 

2 

79.116 

79.118 

82.710 

82.711 

0 

5 

25.158 

25.158 

13.420 

13.420 

15 

5 

25.760 

25  .  764 

15.388 

15.384 

30 

5 

31.399 

31.404 

25.643 

25.644 

45 

5 

44.872 

44.872 

43.711 

43.709 

60 

5 

61.192 

61.189 

62.941 

62.936 

75 

5 

73.902 

73.900 

77.264 

77.262 

90 

5 

78.644 

78.646 

82.534 

82.534 
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THE  LOAD  CORRECTION  PRINCIPLE 


Taft  H.  Broome,  Jr.  Sc.D. 


ABSTRACT 


A  methodology  for  generating  a  geometrically  general  class  of  contlnua  to 
represent  large  lattices  is  presented.  This  methodology  is  characterized  by 
a  mapping  of  the  geometry  and  material  properties  of  the  lattice  into  general 
contlnua  whose  material  properties  can  be  prescribed  independently  of  the  lat¬ 
tice's  architecture,  and  whose  loading  consists  of  the  lattice's  loading 
together  with  an  intrinsic  load  called  the  "load  correction."  This  load 
correction  is  formulated  in  a  manner  that  is  consistent  with  the  constraint 
that  the  continuum  and  the  lattice  exhibit  similar  kinematic  responses  to 
their  respective  loadings.  Knowledge  about  the  lattice  is  obtained  from  a 
Rltz  analysis  of  the  continuum  wherein  efficiency  derives  from  the  repeating 
nature  of  the  lattice's  architecture  and  the  small  set  of  basis  functions 
necessary  to  describe  the  global  response  of  the  continuum. 
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1 .  Introduction 

The  purpose  of  this  study  is  to  assess  the  utility  of  the  load  correction 
principle1  by  which  the  mechanics  of  large  discrete  structural  systems  can  be 
mapped  into  the  mechanics  of  continuous  media.  The  aim  is  to  contribute  new 
experiences  toward  the  production  of  a  comprehensive  theory  that  can  be  imple¬ 
mented  economically  for  the  analysis  of  large  complex  lattice  systems. 

Since  lattices  will  likely  provide  structural  support  for  space  stations, 
large  orbiting  antennas,  solar  power  satellites,  etc.,  then  miles  of  cosmic 
space  may  be  occupied  by  individual  periodic  assemblages  of  structural  ele¬ 
ments  possessing  identical  architectures  (see  Figure  1).  Lattices  of  this 
sort  possess  enormous  numbers  of  degrees  of  freedom  (>104d.o.f).  Such  num¬ 
bers  are  sufficient  to  eliminate  the  finite  element  method  (FEM)  from  conten¬ 
tion  as  an  economical  tool  for  the  global  analysis  of  large  space  structures 
(LSS).  Thus,  more  economical  analysis  methods,  such  as  "continuum  modeling" 
methods,  offer  the  promise  of  economical  alternatives  to  global  application  of 
the  FEM  to  LSS. 

One  class  of  continuum  modeling  methods  focuses  upon  similarities  -  where 
they  exist  -  between  the  material  domains  of  lattices  and  simple  continua  such 
as  rods,  beams,  shafts,  plates,  and  shells  .  Material  properties  for  a  con¬ 
tinuum  are  sought  that  enable  it  and  the  lattice  it  attempts  to  model,  to 
exhibit  similar  kinematic  responses  to  similar  loadings  while  being 
constrained  by  similar  displacement  boundary  conditions.  Thus,  the  large  set 
of  algebraic  equations  associated  with  FEM  models  of  LSS  are  avoided  in  favor 
of  small  sets  of  differential  equations  to  which  closed-form  solutions  exist, 
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or  conventional  solution  techniques  are  applicable.  The  load  correction 
method  (LCM),  which  is  a  means  of  implementing  the  load  correction  principle, 
is  a  straight  forward  general  continuum  modeling  method  which  is  adaptable 
to  both  simple  and  complex  domains. 

The  LCM  permits  a  priori  specifications  of  the  material  properties  of  the 
continuum  model,  and  seeks  a  new  loading  for  it  that  is  consistent  with  the 
constraint  that  the  model  and  the  lattice  exhibit  similar  kinematic  responses 

3 

(  including  observance  of  similar  kinematic  displacement  boundary  conditions) . 
As  we  shall  see,  this  new  loading  consists  of  the  loading  applied  to  the  lat¬ 
tice  and  a  displacement-dependent  load  called  the  ’’load  correction.**  Because 
closed-form  solutions  are  not  yet  readily  apparent,  a  Ritz  analysis  is  used  to 
transform  the  large  set  of  algebraic  equations  produced  by  the  FEM  into  a 
smaller,  more  manageable  set.  Also,  the  inner-product  nature  of  a  Ritz  ana¬ 
lysis  enables  a  Dirac-delta  technique  to  take  advaptage  of  the  repeating 
architecture  of  the  LSS  to  arrive  at  this  small  set  of  algebraic  equations 
economically. 

The  objectives  of  this  study  are  as  follows: 

(a)  Introduce  a  nomenclature  that  facilitates  the  uncoupling  of  globally 
dependent  parameters  from  locally,  or  what  will  be  defined  as  ” cluster"  depen¬ 
dent  parameters; 

(b)  Formulate  the  given  lattice  problem; 

(c)  Map  the  given  lattice  problem  into  a  continuum  model  according  to  the 
load  correction  principle; 

(d)  Solve  the  mod*’ ; 

(e)  Illustrate  how  the  globally  dependent  parameters  can  be  evaluated 
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economically;  and 

(£)  Demons Crate  che  utility  of  the  LCM  via  sample  problems. 

The  utility  of  the  LCM  will  be  assessed  in  terms  of  comparisons  with  the  FEM 
and  an  LCM/FEM  hybrid  approach;  and  in  terms  of  the  degree  to  which  the  appli¬ 
cations  spectrum  of  the  LCM  extends  beyond  those  of  other  continuum  modeling 
methods . 

The  scope  of  this  study  is  limited  to  the  static  analysis  of  planar 
lattices.  This  is  the  narrowest  scope  within  which  the  objectives  of  this 
study  can  be  realized,  and  it  is  adequate  as  an  introduction  to  the  LCM. 
Moreover,  even  though  modal  analysis  of  LSS  is  of  immediate  interest  to  dyna¬ 
mic  control  efforts,  the  role  of  the  LCM  in  these  efforts  will  likely  take 
che  form  of  a  static  tool  for  generating  generalized  or  reduced-order  flexibi¬ 
lity  matrices  for  LSS.  Thus,  this  study  is  the  first  stage  in  the  development 
of  a  comprehensive  LCM  theory. 

The  remainder  of  this  study  is  organized  into  seven  sections  and  three 
appendices.  The  following  six  sections  address  the  six  objectives  respec¬ 
tively,  and  the  seventh  summarizes  them  with  conclusions  about  the  utility  of 
the  LCM.  The  appendices  contain  a  list  of  the  LCM's  nomenclature,  definitions 
of  some  analysis  parameters,  and  a  list  of  references. 

2.  The  Cluster  Concept 

The  given  lattice  structure  will  be  modeled  by  a  continuum  whose  rnterial 
domain  spans  a  configuration  of  points  that  is  always  identical  to  the  con¬ 
figuration  of  the  lattice's  joints  (for  example,  see  Figures  2a  and  b) .  These 
points  in  the  continuum  are  called  the  "images"  of  the  joints  in  the  lattice. 
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The  index  "y"  identifies  each  joint  in  the  lattice  as  well  as  its 
corresponding  image  in  the  continuum.  Augmenting  this  global  numbering  system 
is  a  nomenclature  that  associates  each  joint  and  its  image  with  locally 
repeating  architectures  called  "clusters.*’ 

A  lattice  cluster  is  defined  as  a  subdomain  of  the  given  lattice  that 
facilitates  the  construction  of  a  stiffness  equation  relating  the  external 
forces  at  a  joint  to  the  translational  displacements  of  this  joint  and  of 
other  neighboring  joints.  This  locally  dependent  equation  will  be  called  the 
"cluster  stiffness"  equation.  From  the  examples  of  lattice  clusters  shown  in 
Figures  2c  and  d,  it  is  clear  that  the  clusters  are  themselves  assemblages  of 
repeating  elements  (or  cells)  that  have  been  the  focus  of  analysis  in  other 
continuum  modeling  approaches. 

Similarly,  the  images  of  the  joints  defining  a  lattice  cluster  will  define 
its  associated  continuum  cluster. 

A  set  of  identical  lattice  clusters,  each  being  associated  with  identical 
continuum  clusters,  are  identified  collectively  by  the  index  "t".  In  Figure 
2c,  for  example,  the  dusters  bounded  by  global  sets  {y  -  1,2,3,51,52,53}, 

{y  *  3,4,5,53,54,55},  etc,  are  identified  as  clusters  of  type  t  ■  1.  Sets 
(y  •  2,3,4,52,53,54},  {y  ■  4,5,6,54,55,56},  etc.,  are  identified  as  clusters 
of  type  t  -  2. 

With  respect  to  each  cluster  type,  the  joints  and  images  at  which  the  for¬ 
ces  in  the  duster  stiffness  equations  are  applied  are  called  "principal" 
joints  and  images,  respectively,  in  the  cluster,  and  are  identified  by  the 
index  "e."  Moreover,  the  joints  and  images  at  which  the  displacements  in 
these  equations  are  required  are  identified  by  the  index  "g."  The  index  "a" 
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is  used  Co  distinguish  two  or  more  principal  joints  or  images  sharing  the  same 
(x,e)  designation,  but  different  y  designations.  We  define  a  so  that  there  is 
a  one-to-one  relationship  between  each  y  and  each  3-tuple  (x, e, a). 
Symbolically, 


y  -  §(x,e,a) 


(x,e,a)  ■  9  *(y) 


2.1 


where:  9  is  an  operator;  and  9  ^  is  its  unique  inverse  (see  Table  1). 

The  4-tuple  (x,e»a,B)  identifies  a  unique  point  y  in  a  global  domain 
according  to  the  rule 

Y  *  C<Tte,a,8)  ,  (x,s,a, 0)  •  6  *(y)  2*2 

where:  6  is  an  operator;  and  its  inverse  G  *  exists,  but  is  not  necessarily 
unique.  We  define  G  so  that  the  spadal  relationship-  shared  by  points 
f(T,e,a)  and  G(T,e,a,0)  in  the  global  domain  is  Identical  to  the  spatial 
relationship  shared  by  points  identified  by  the  pairs (t,e)  and  ( x,  0)  in  the 
clusters,  respectively  (see  Table  2). 

Figures  2c  and  d,  together  with  Tables  1  and  2,  define  a  cluster 
"scenario"  by  which  locally  dependent  stiffness  equations  can  be  constructed 
at  each  joint  in  the  lattice,  and  stiffness  equations  of  similar  forms  can  be 
constructed  at  the  images  of  these  joints  in  the  continuum.  Other  scenarios 
are  possible  -  some  affording  increased  precision  in  the  cluster  stiffness 
equations,  but  at  the  expense  of  increased  computational  requirements. 

In  the  analysis  which  follows,  the  Greek  indices  x,  e,  a  and  g  are  used  as 
superscripts  to  distinguish  discrete  parameters  from  their  continuous  forms 
(see  Table  3  ).  Roman  subscripted  indexes  m,  n,  i  and  j  denote  components  of 
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TABLE  1 :  The  9 -transformation  for  Figure  2. 


1 

1 

2 

3 

24 

1 

KUDU 


2  2  2  2  2  2  2  2  2  2  2  2 


1  1  2  3  -  24  1  1  2  3  -  24  1  2  3  -  24  1  1  2  3  -  24 


IDQE 


1  2  4  6  -  48  51  52  54  56  -  98  3  5 


E 

9 


ro 

2 

2 

i 

1 

1 

i 

2 

3 

49  50  53  55  57  -  99  100 


TABLE  2:  A  sample  of  the  ©-transformation  for  Figure  2. 


III 

1 

1 

2 

1  2  2  2  2  2  2  2  2  2  2  2  2 


G 

G 

IDEGGI 


123456123456 


1  1  2  3  51  52  53  3  4  5  53  54  55 


I 

I 

Hi 

IEEDI 

Hi 


.  1 

1 

2 

2 

2 

2 

5 

6 

TABLE  3;  Continuous  and  discrete  parameter  forms. 


Parameter 

Type 

Continuous 

Forms 

Discrete 

Forms 

Cartesian 

coordinates 

X 

n 

xY 

n 

Tea 

X 

n 

*tcaB 

n 

Position 

Vectors 

Hi 

r  ■  Z  x  e 
n  n  n 

-•■Tea  _  Tea  ■* 

r  »  Z  x  e 

n  n  n 

jTeaB  .  xTeaB  S 
n  n  n 

Displacement 
Basis  Functions 

\  ■  v*> 

Tea  ,+Tea. 

Sn  "  Vr  ) 

Tea8  .-►TeaS. 

"  qm(r  } 

Assumed  Dis¬ 
placements 

%  ■  un<?) 

un  -  un<r) 

Tea  „  (rTea) 

n  n 

uTEa6  -  u  (Vta&) 
n  n 

Displacement 

Vectors 

-*■  -*■ 
u  *  Z  u  e 
n  n  n 

■*Y  ^  Y  -*■ 

u  *  l  u  e_ 

n  n  n  n 

-►Tea  _  Tea  -*■ 

u  =  Z  u  e 

n  n  n  n 

-►TeaS  _  TeaS  -► 

u  *  i  u  e 

n  n  n  n 

40 


8 


these  parameters;  and  m' ,  n' ,  i'  and  j'  are  alternates  of  their  unprimed  coun¬ 
terparts.  The  value  of  this  nomenclature  is  that  it  facilitates  the 
uncoupling  of  global  from  local  dependencies  since  parameters  that  are  depen¬ 
dent  upon  a  and/or  y  are  globally  dependent;  and  those  that  are  Independent  of 
a  and/or  y,  but  dependent  upon  t,  e  and/or  g,  are  locally  or  cluster  depen¬ 
dent. 


3.  The  Given  Lattice  Problem 

The  lattice  cluster  stiffness  equations  are  given  in  the  form 


£  £ 
n*  8 


rreB 

nn* 


_teoB 


(t,e,a)  -  9_1(y) 


3.1 


where:  the  overbar  refers  to  the  lattice;  Is  the  component  of  the  external 

force  acting  or*  the  lattice  at  joint  y  in  the  direction  of  the  global  unit 

base  vector  e  «  k™?  are  the  cluster  stiffness  coefficients,  or  simply  the 
n  nxi 

"stiffnesses;*'  and  u^f°®  are  the  displacements  of  the  joints  in  the  cluster. 

In  the  case  that  the  lattice  is  a  truss  (i.e.,  all  joints  are  pinned), 
equation  3.1  is  the  exact  form  of  the  static  equilibrium  equation. 

If  the  lattice  is  a  frame  (i.e.,  all  joints  are  rigid)  or  a  hybrid  of 
pinned,  rigid  and/or  flexible  joints,  the  equilibrium  equations  can  be  reduced 
to  the  form  of  equation  3.1  by  (a)  invoking  St.  Venant's  principle  to  neglect 
the  Internal  moments  at  joints  6  in  the  cluster  which  are  remote  from  the 
principal  joint  "  e  ";  and  (b)  using  a  standard  condensation  procedure  to  elimi¬ 
nate  the  external  moment  equations. 

We  shall  assume  that  for  each  pair  (n,y)  where  is  unknown,  the  displa¬ 
cements  are  prescribed  ns 
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F-F  3.2 

n  n 

where:  are  known.  The  problem  Is  to  find  that  satisfies  equation  3.1 

wherever  F  is  known, 
n 

4.  The  Continuum  Model 

A  continuum  will  be  said  to  be  "topologically'*  equivalent  to  a  lattice  if 
the  undeformed  material  domain  of  the  continuum  spans  a  set  of  material 
points,  called  "images,"  whose  geometrical  configuration  is  identical  to  the 
undeformed  configuration  of  joints  in  the  lattice.  A  continuum  that  is  topo¬ 
logically  equivalent  to  a  lattice  will  also  be  said  to  be  "kinematically" 
equivalent  to  it  if  all  deformed  configurations  of  the  lattice's  joints  and 
their  images  in  the  continuum  are  Identical.  As  we  shall  see,  the  choices  of 
material  domain  and  displacment  constraints  at  non-image  points  of  the  con¬ 
tinuum  will  be  motivated  by  factors  of  simplicity. 

A  continuum  model  for  a  lattice  is  defined  as  a  continuum  that  is  kinema¬ 
tically  equivalent  to  the  lattice.  The  formulation  of  the  model  begins  with  a 
choice  of  material  domain  that  is  topologically  equivalent  to  the  lattice,  and 
a  choice  of  material  or  constitutive  properties  that  may  be  specified  indepen¬ 
dently  of  the  lattice.  This  choice  will  also  be  motivated  by  factors  of 
simplicity. 

Static  equilibrium  on  the  surface  S  of  the  continuum  model  is  satisfied  by 

T  (r)  -  d  (u)  4.1a 

n  n 

where:  T  is  the  component  of  the  surface  traction  in  the  e  direction;  and 
n  n 

dn  is  a  scalar  differential  operator.  Static  equilibrium  in  the  volume  V 
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enclosed  by  S  is  satisfied  by 


B  (r)  -  D (u)  4.1b 

n  n 

where:  Bq  is  the  component  of  the  body  force  in  the  eQ  direction;  and  Dq  is  a 
scalar  differential  operator. 

Let  SY  be  any  subdomain  sts  containing  only  one  image  yeSY.  Let 
6Wy(Tq)  be  the  variation  of  worlc  done  on  the  model  by  all  the  surface  trac¬ 
tions  Tq  acting  on  SY.  Then 


«W 


(Tq) 


;/ 


6u_  dS 


Similarly, 


6WY(Bn) 


XJJ 

VY 


4.2a 


4.2b 


where:  6Wy(B  )  and  VY  are  work  and  volume  domain  analogues  of  5WY(T  )  and 
n  n 

SY  respectively.  These  conditions  must  clearly  be  satisfied  no  matter  how  the 
loads  on  the  continuum  are  applied. 

We  shall  stipulate  that  the  loads  shall  be  applied  so  as  to  satisfy 


6Wy(T  )  -  Fy  5uY  ...  -V-rcS  ,  and  4.3a 

n  n  n  9 

6Wy(Bq)  -  FY  5uY  ...  -V-ycV  4.3b 
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for  all  pairs  (n,y)  where  are  prescribed,  and  where: 


-or  .  -  -  tTe8  reoB  fTY  _  _  r-reB  teaB 

Fn  EZ  knn*  V  +  Fn  Z,  \  knn‘  V  » 

n'B  n  B 


I,  Z  k «8  uTfaB  6uy 
n  g  nn  n  n 


II 


d_(u)  5u„  dS  ...  yeS 
n  n 


III 


Dq(u )  5uq  dV  ...  y  c  V  . 


The  coefficients  can  be  obtained  from  finite  difference  approximations  of 

the  right-hand  sides  of  equations  4.2a  and  b.  Below,  these  stiffness  coef¬ 
ficients  will  be  obtained  from  use  of  the  cluster  as  a  device  for  imple¬ 
menting  the  finite  element  method  (FEM)  on  a  local  scale. 

In  either  case,  substitution  of  equation  4.4a  into  equations  4.3a  and  b, 
and  appeal  to  equations  4.2a  and  b,  give 


II 


d  (u)<5u  dS 
n  n 


<E  E  kT‘fuTf«B  +  71  -  E  Ik!‘f 

,  _  nn  n  n  ,  _  nn  n  n 

n'  B  n'  B 


where  ycS,  and 


III 


l  (u)6u  dV  -  (I  E  ♦  T  -  E  E  uTf»8)Sar 

n  n  ,  nn  n  n  ,  nn  n  n 

n'  8  n'  B 


where  ^«V.  Then,  in  view  of  equation  4.4b,  we  can  deduce  from  equations  4.5a 
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and  b  Chat 


-y_v  ufe8 

“  t  nn'  V 
n’  S 


for  all  pairs  (n,y)  where  SuY  is  arbitrary  (i.e.,  where  are  prescribed), 


By  further  stipulating  that 


uY  -  F 
n  n 


for  all  pairs  (n,y)  where  tT[  are  prescribed  (l.e.,  where  6uY  ■  0),  we  can 


conclude  that 


uY  -  u^  ...  ^Kn,y) 


Thus,  the  continuum  Is  kinematically  equivalent  to  the  lattice. 

As  long  as  the  axiom  of  continuity  is  not  violated  (i.e.,  the  Jacobian 
1*0),  kinematic  equivalence  Is  preserved  regardless  of  any  additional  displa¬ 
cement  constraints  imposed  at  non-image  points  In  the  continuum.  Here  again 
simplicity  will  be  seen  as  the  motivating  factor  for  imposing  displacement 
constraints  on  the  continuum  at  non-image  points. 

Let  the  "load  correction"  be  defined  as  the  force 


r  r  AltTEB  ,.Tect0 

Z,  l  ^nn*  V 

n'  0 


where 
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Ak 


T£0 

nn' 


, -reS  _  -r teB 
on'  nn' 


4.9b 


Then 


fy  -  F  + 

n  n 


AF ' 


4.10 


In  summary: 


A  continuum  that  is  topologically  equivalent  to  a  lattice 
is  also  kinematically  equivalent  to  it  if  (1)  the  displace¬ 
ment  prescriptions  imposed  on  the  lattice  are  imposed  on 
the  continuum  at  their  corresponding  images;  and  (2)  the 
load  on  the  continuum  is  the  sum  of  the  lattice  forces  and 
the  load  corrections . 


This  principle,  called  the  "load  correction  principle,”  Is  the  strategy  by 
which  the  given  lattice  problem  was  mapped  into  a  continuum  model. 

3.  Solution  of  the  Model 

An  approximate  global  solution  of  the  continuum  model  will  be  obtained 
using  a  Rltz  approach.  The  basis  functions  for  the  assumed  displacement  field 
are 

«■  “  <m(V  >  5*2 

u  TGg 

where:  the  range  of  m  is  finite;  and  the  various  discrete  forms  q^, 
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and  q^Ect®  are  given  in  Table  3.  To  simplify  the  analysis,  two  special  con¬ 
ditions  will  be  imposed  on  this  set  of  basis  functions: 


TfoS  _ 


9V3xn 


Tea 

V 


D 

nmm 


'^m' 


f 


and 


5.2a 

5.2b 


T  c  6 

where:  J _ ,  and  D _ .  are  known  constants. 

This  set  if  basis  functions  will  generate  the  assumed  displacement  field 
in  the  fora 

un  S  bna^m  +  Z  f  anmj  V  Xj  5*3 

m  ■  j  J  J 

where:  u  is  defined  in  Table  3;  and  b  and  a  .  are  coefficients  that  are 
n  nm  nmj 

specified  so  that  the  displacement  boundary  conditions  are  satisfied  for  all 
choices  of  the  varying  parameters  Then,  the  various  discrete  forms  of 

un  shown  in  Table  3  can  be  expressed  in  terms  of  \ ^ .  Furthermore,  the 
variation 


5un 


l  5 

m  J 


5.4 


and  its  various  discrete  forms  follow  as  suggested  by  Table  3. 

Using  these  equations,  the  forces  applied  to  the  continuum  can  be  written 
as 
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fy  -  r  ♦  r  b**Tve6  +  z  z  a**je  <ca  ^ 

n  n  ns  in  ,  naj  to  j 

a  a  j 


where:  the  constants  bQffl  and  anoj  are  given  in  Appendix  2.  The  variation 

of  work  on  the  model, then,  is 


5W  -  Z  Fr  6uy 
n  n 
Y 


which,  after  substitutions  from,  equations  5.5  and  5.3,  becomes 


«w  -  *  ♦j(b)  ♦  )5XJ 


where:  •j^,  and  *j*«  are  defined  in  Appendix  2. 

Using  engineering  notation,  the  components  of  strain  in  the  continuum  are 


Ei  ’  S  Z,Ginn»  7“ 
n  n'  3x  , 

n 


in  the  linear  case  where:  ^nn>  ere  known  constants.  Differentiation  of 
equation  5.3,  and  substitution  of  equation  5.2b  into  -he  result  yield 


ei  ■ 1  <s£>+  f 

o  j 


where:  G^and  are  given  in  Appendix  2. 

If  we  assume  that  the  model  obeys  Hooke's  law,  the  stress  can  be  written 
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aL  “  J,  Cii,ei' 


5.9 


where:  C^,  are  prescribed  consdcuclve  constants.  Then,  the  variation  of 
internal  energy 


5P  -  l 


dV 


5.10 


can  be  written  as 


6P 


’  ] ( J.  fjjV5xj 


5.11 


after  appropriate  substitutions  of  equations  5.8  and  5.9  into  equation  5.10, 
(b)  (&) 

and  where:  Y  j  and  Yj ^ ,  are  given  in  Appendix  2. 

Since  SP  «  SW,  and  since  6Xj  are  arbitrary,  the  solution  for  Xj  proceeds 
froa  (a)  setting  equation  5.6b equal  to  equation  5.11  to  get 


,  ,.<*>  .  .(*)■,,  .  .0>)  .  ,0>>  .  .(F)  . 

),  CTJJ'  ’j  *1  *3  • 


5.12 


and  (b)  implementing  the  standard  Gaussian  elimination  procedure.  Back 

substitution  of  the  resulting  Xj  into  equation  5.3  gives  the  continuous 

displacement  function  for  the  model,  and  the  discrete  form  u  Y  is  the  solution 

m 

of  the  lattice. 


6.  Computational  Economy 


The  effort  required  to  obtain  an  FEM  solution  to  the  given  lattice  problem 
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can  be  measured  roughly  In  terms  of  the  cpu  (and  memory)  required  to  solve  a 
narrowly  banded  set  of  algebraic  equations  of  order  TN  where:  1  <  y  <  n  and 
1  <  n  <  N.  The  effort  required  to  obtain  and  LCM  solution  to  the  same  problem 
can  be  measured  roughly  in  terms  of  the  cpu  required  to  (a)  solve  a  dense  set 
of  algebraic  equations  of  order  J,  where  1  <  j  <  J  (see  equation  5.12);  and  to 
(b)  evaluate  the  parameters  listed  in  the  Appendix.  As  we  shall  see,  J  «  fK 
generally,  and  the  evaluation  of  these  parameters  can  be  made  with  an  effort 
that  is  generally  independent  of  the  ranges  of  the  global  indexes  a  and  y. 

From  Appendix  2,  it  is  clear  that  if  the  computational  effort  required  to 
evaluate  the  parameters  4^  and  0T£ ,  is  independent  of  A  (where  1  <  a  <  A)  and 
F,  then  the  effort  required  to  evaluate  the  remaining  parameters  is  also  inde¬ 
pendent  of  A  and  r.  To  see  what  effort  is  actually  required  to  evaluate  these 
parameters,  consider  the  equation 

-  I  ?!  ql  .  6.1 

Tim  Q  Til 

y 

An  important  role  of  the  LCM  in  the  analysis  of  large  space  structures  is 
to  produce  a  generalized  or  reduced-order  flexibility  matrix  for  the  lattice 
so  that  the  frequencies  and  mode  shapes  can  be  deduced  from 

[I  -  oi2fM]fu}  -  0  6.2 

where:  I  is  an  identity  matrix;  u>  is  a  natural  frequency;  f  is  the  flexibi¬ 
lity  matrix,  which  is  the  inverse  of  the  stiffness  matrix;  M  is  the  genera¬ 
lized  mass  matrix;  and  y  is  the  mode  shape  vector.  The  flexibility  matrix 
will  be  generated  column-wise  from  sequential  applications  of  a  single  unit 
force  to  the  model.  Thus,  each  use  of  equation  6.1  will  usually  involve  only 
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one  non-zero  value  of  F.  This  is  to  say  that 

U 

*nm  "  F”o  q®o  ^n0n^m0m  •••  Q0  8um  on  n0  and  °0  *  6.3 

where:  the  Indexes  n  and  m  are  fixed  values;  and  6-  n  and  6™  m  are 

O  O  uO  1  “*0*“ 

Kronecker  deltas. 

If  gravity  gradients,  solar  pressure  or  other  distributed  loads  are 
applied  to  the  lattice,  then  equation  6.1  can  be  computed  economically  if  we 
can  write 


T  -  Z  c  qT 
n  nmMm 


6.4 


where:  c  are  known  constants.  Substitution  of  this  form  into  equation  6.1 
am 

yields 


♦o.  •  C,  =n.'  <5  • 

m  Y 


6.5 


The  term  in  parentheses  can  be  evaluated  efficiently  using  a  Dlrac-delta 


technique. 


For  example,  if 


q 


Y 

2 


q2 


(xj ) 


Y>2 


6.6 


then 


I 

Y 


Y  Y 
q2  q2 


£ 

Y 


6.7 


But,  using  the  Dirac-delta  function,  we  can  write 
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(X  ) 


-  J  52  6(z-xj) d£ 


6.8 


so  chat 


+x[+0 

1  q2  q2  *  I  ?2  tz  s<5-*J)]de  6.9 

Y  J  i  Y 


where  the  term  in  brackets  is  an  integrable  function.  Integrating  equation 
6.9  by  parts  gives 


r  „Y  _Y 
^2  ^2 


U2z, 


sz2  +  Z3]  1 


6.10 


where:  -Jl  «(§-x"J)dC;  and  Zfc+1  -  J*  Z^dS,  k-1,2.  The  integrals  Zj,  Z2, 

can  be  integrated  exactly  -  by  hand.  Back  substitution  of  these  values 
into  equation  6.10,  and  then  into  equation  6.5  results  in  a  symbolic  eva¬ 
luation  of  Thus,  the  cpu  required  to  evaluate  #  is  independent  of  r, 

and  the  memory  requirements  are  virtually  miniscule. 

In  a  like  manner,  the  parameter 


TC 


a 


tea  tea 


V 


6.11 


is  evaluated  using  the  above  Dirac-delta  technique,  thereby  implementing  a 
globally  independent  effort.  Thus,  the  LCM  can  be  implemented  economically. 
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7 .  Sample  Problems 

A  selection  of  plane  lattice  structures  is  discussed  below  as  a  means  of 
demonstrating  the  versatility  of  the  LCM.  In  each  case,  the  continuum  model 
is  an  isotropic,  plane  stress  solid  whose  material  domain  follows  the  outline 
of  the  lattice.  Where  the  displacement  boundary  conditions  for  the  lattice 
are  discrete  hinge  supports,  the  model  is  constrained  to  obey  hinge  support 
conditions  along  line  segments  joining  the  images  of  these  hinges.  The 
assumed  displacement  field  is  generated  by  power  basis  functions 

{q^:  1,  Xj,  Xj,  (xj)  ,  »^x2^  »  x2>  •••  }•  ^.1 

The  LCM  was  implemented  via  the  LCM21  computer  code  on  a  Vax  11/750  ,  and  the 
FEM  was  implemented  via  STRUDL  on  an  IBM  30335. 

Problem  #1.  TWo  loading  cases  are  considered  for  the  cantilever  truss 
shown  in  Figure  3  .  The  truss  is  a  single-braced,  one-way,  single-story 
structure  that  is  well  known  to  behave  much  like  an  orthotropic  beam.  The 
constitutive  properties,  however,  will  be  modeled  as  isotropic  in  the 
following  LCM  model. 

In  the  boom's  case  ,  the  LCM  models  the  horizontal  displacement  field  with 
a  linear  displacement  field.  This  model  is  sufficient  no  matter  what  the 
length  of  the  boom  may  be.  Thus,  the  cpu  required  by  the  LCM  is  Independent 
of  the  global  d.o.f.  of  the  boom  (see  Problem  #2  below). 

In  the  bending  case,  the  LCM  models  the  horizontal  displacement  field  with 
a  parabola,  and  the  vertical  with  cubic  parabola.  Comparison  with  the  FEM  is 
adequate,  and  cpu  remains  independent  of  the  global  d.o.f.  of  the  system. 

Problem  #2 .  The  boom-like  loading  case  is  considered  for  the  cantilever 
lattice  shown  in  Figure  4  .  Two  cases  are  considered:  truss  (or  pinned 
jointed);  and  frame  (or  rigid  jointed). 
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A  linear  displacement  field  is  used  to  model  the  horizontal  displacements, 
and  various  orders  of  parabolic  fields  are  considered  for  the  vertical  displa¬ 
cements.  A  cluster  of  two  panels  was  sufficient  for  the  truss,  and  five 
panels  were  used  in  the  case  of  the  frame.  While  increasingly  larger  displa¬ 
cement  fields  for  v  are  required  to  capture  the  edge  effects  in  increasingly 
longer  booms,  a  linear  displacement  field  remains  sufficient  for  u  whatever 
the  overall  length  of  the  boom*  Because  the  d.o.f.  of  these  displacement 
fields  are  small,  the  LCM  cpu  is  virtually  Independent  of  the  length  (i.e., 

#  bays),  while  the  FEM  cpu  increases  geometrically. 

Problem  #3.  The  boom-like  complex  case  is  considered  for  the  doubly- 
braced  hybrid  of  pinned  and  rigid  joints  shown  in  Figure  5  .  The  horizontal 
displacement  field  is  linear,  and  the  vertical  is  of  10  d.o.f. 

The  LCM  and  FEM  results  for  the  horizontal  displacements  are  almost 
indistinguishable.  The  vertical  displacements,  however,  compare  favorably  in 
magnitude,  sign,  and  overall  trends.  Edge  effects  are  captured  to  the  degree 
allowed  by  the  assumed  displacement  field  (Figure  S  ). 

Problem  #4.  Shown  in  Figures  6-8  are  two-story  cantilever  lattices 6. 

Three  local  architectural  patterns  are  considered:  doubly-braced,  diamond 
pattern  of  alternating  diagonals,  and  W-pattern  of  alternating  diagonals. 
Comparisons  of  FEM  with  LCM  results  are  shown,  and  these  comparisons  generally 
follow  those  observed  in  the  single-story  cases:  longitudinal  displacements 
for  booms  approximated  well  in  the  linear  case;  bending  results  were  good; 
etc.  There  were,  however,  some  exceptions. 

The  bending  cases,  for  relatively  short  beams,  required  higher-order 
approximations  for  the  displacement  fields,  and  cpu  for  the  short  beams  cora- 
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FIGURE  7:  Two-story,  d iamond -braced  cantilever  trusses. 
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pared  with  FEM  cpu.  However,  FEM  cpu  was  noc  competitive  with  LCM  cpu  for 
longer  beams. 

The  vertical  displacements  shown  in  Figure  7  were  not  approximated  well. 
Although  horizontal  displacements  were  good  in  this  case,  as  much  as  12Z  error 
in  these  displacements  was  observed  in  the  single-story  case.  However,  at  the 
expense  of  non-competitive  cpu,  the  LCM  accurately  predicted  the  displacement 
behavior  of  short  booms  assembled  from  this  pattern  in  both  the  single-story 
and  the  two-story  cases7.  The  failure  of  the  LCM  in  the  longer  cases  is, 
therefore,  attributable  to  the  choice  of  assumed  basis  functions  of  the  power 
type  where  the  actual  vertical  displacement  field  varies  periodically  in  the 
same  period  as  the  local  repeating  architectures  (see  Figure  7  ).  A  Fourier 
basis  set  may  well  resolve  this  problem. 

Problem  #5.  The  flat,  ln-plane  deforming  plates  shown  in  Figure  9  are 
considered.  Singly-braced,  alternating  diagonal  trusses  are  considered,  and  a 
hybrid  LCM/FEM  analysis  Is  compared  with  a  conventional  global  FEM  analysis. 

The  problem  is  to  compute  the  vertical  displacement  under  the  load  applied 
to  the  5-bay  x  5-bay  plate.  Two  convergent  cases  are  considered.  First  ,  the 
FEM  is  applied  to  Increasingly  larger  subdomains  of  the  plate  assuming  the 
internal  joints  are  pinned. 

Second,  the  LCM  is  used  to  commute  the  displacements  of  the  5-bay  x  5-bay 
plate.  These  results  are  then  used  as  boundary  conditions  at  the  internal 
joints  where  the  FEM.  is  again  applied. 

The  second  case  converges  much  faster  than  the  first,  producing  an  18Z 
improvement  in  system  order  at  the  15  d.o.f.  subdomain.  Thus,  the  hybrid 
LCM/FEM  is  more  economical  than  the  conventional  global  FEM  for  this  case. 
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FIGURE  9:  In-plane-deforming  plate-like  lattice. 
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8.  Summary  and  Conclusions 

A  nomenclature  was  introduced  in  section  2  that  facilitated  the  uncoupling 
of  duster  dependent  parameters  from  globally  dependent  ones.  The  given  lat- 
tice  problem  was  then  formulated  in  terms  of  cluster  coordinates.  A  continuum 
model  was  formed  as  a  continuum  that  is  topologically  and  kinematically 
equivalent  to  the  lattice,  but  subject  to  a  displacement-dependent  loading 
called  the  “load  correction."  Using  this  load  correction,  the  continuum  model 
was  described  In  terms  of  cluster  coordinates. 

The  model  was  solved  using  a  Ritz  procedure  which  introduced  globally 
dependent  parameters.  However,  when  considering  lattices  with  repeating 
cluster  architectures  -  as  will  typify  LSS  -  these  parameters  were  reduced  to 
cluster  dependencies  using  a  Dirac-delta  technique  for  computing  long  sum¬ 
mations.  Sample  problems  were  presented  which  demonstrated  that  the  range  of 
lattice  types  to  which  the  LCM  is  adaptable  exceeds  the  present  scope  of  com¬ 
peting  continuum  modeling  methods.  While  these  methods  lead  to  closed-form 
solutions,  the  LCM  leads  to  numerical  models. 

Some  Important  features  of  the  LCM  are: 
e  General  applicability  to  pin-jointed,  rigid-jointed,  and  hybrids  of 
pinned-  and  rig id- jointed  lattices; 

•  Mlcropolar  models  are  not  necessary  for  rigid-jointed  lattices  as  the 
cluster  is  a  means  of  invoking  St.  Venant's  principle  locally  so  that  the 
rotational  degrees  of  freedom  can  be  eliminated  via  conventional  static 
condensation  processes; 

•  The  matrices  comprising  the  LCM  are  either  the  familiar  FEM  stiffness 
matrices  whose  dimensions  are  of  the  cluster  d.o.f.,  or  are  otherwise 
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easily  generated; 

•  The  LCM  can  be  applied  to  a  vide  range  of  lattices  possessing  different 
local  and  global  architectures  without  adjusting  the  methodology  to  accom¬ 
modate  them; 

•  The  LCM  accommodates  a  choice  of  assumed  displacement  fields.  Although 
present  considerations  Included  power  basis  functions,  Fourier  basis 
function  may  be  desirable  In  cases  where  the  displacements  of  the  truss 
vary  periodically  In  the  same  period  as  the  clusters; 

•  The  LCM  has  the  ability  to  capture  edge  effects; 

e  The  global  response  of  the  LCM  can  be  captured  at  a  cpu  that  Is  substan¬ 
tially  less  than  that  of  FEM  for  large  d.o.f.  structures.  In  most  cases 
observed,  the  LCM  cpu  is  virtually  independent  of  the  d.o.f. 

•  Hybrid  LCM/FEM  promises  to  serve  as  a  comprehensive  analysis  tool  for  LSS 
that  is  superior  to  global  Implementation  of  the  FEM. 

Finally,  well-known  generalized  or  lumped  mass  techniques  provide  the 
basis  for  modal  analysis  of  LSS.  The  LCM  or  hybrid  LCM/FEM  can  be  used  as  an 
economical  means  of  obtaining  the  generalized  flexibility  matrices  for  LSS. 
Thus,  a  reduced  basis  dynamic  analysis  technique  for  LSS  proceeds  from  use  of 
the  LCM  as  a  static  tool. 
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APPENDIX  1  :  NOMENCLATURE 
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APPENDIX  2  :  MAIN  PARAMETER  DEFINITIONS 
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ABSTRACT 


This  paper  addresses  computational  Lssues  associated 
with  a  function  D:  I  x  I  — *  R  defined  by  the  relation 


0(N,  I)  -  h 


Z  fV°  f'1  -  \'2 

n-l  J  Xj-0  J  *,-0  J  *i-0 


)  d*, . . .d* 
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*  j  1  *  *1  r  1 

where:  N  and  I  are  positive  definite  integers;  0<h;  whenever 
1  <  N,  we  have  x  .—  x  •  h  for  every  linsN-1;  and  i  is 
Dirac's  operator.  1  n 

INTRODUCTION 


The  motivation  for  this  paper  is  the  appearance  of  D(N,I) 
in  the  LCM. FORTRAN1  algorithm  which  enables  macro¬ 
engineering  analyses  of  structures  exhibiting  local  or  grain- 
size  effects  above  the  micro-scale1  .  In  the  application  of 
the  LCM  to  a  given  mesomechanics1  problem,  we  can  anticipate 
the  integer  I  to  remain  fixed  at  some  number  I  -  100,  while  N 
may  take  on  a  dozen  (or  fewer)  values  spanning  1  to  maxN  » 
10*  .  All  of  the  appearances  of  N  in  the  LCM  algorithm  are 
contained  in  D(N,I).  Therefore,  we  seek  to  studv  computational 
issues  that  are  contingent  to  effective  computer  resources 
management  of  equation- l  when  applying  the  LCM  to 
mesomechanics  problems. 

_  To  sharpen  the  focus  of  this  study,  we  cal!  an  arithmetic 
equivalent  or  ’'exact'*  arithmetic  restatement  of  the  right- 
hand-side  of  equation-,  a  "computable  fomulation  of  .f 

such  a  restatement  can  be  said  to  be  readily  translatable  uti 
FORTRAN  statements.  In  restricting  the  present  discussion  to 
"exact"  formulations  -  rather  than  including  approximations 
we,  thereby,  define  the  first  computational  issue  of  the 
study.  For  example,  a  computable  formulation  produced  bv 
integrating  the  right-hand-siie  of  equation-1  is 


where : 
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The  second  issue  concerns  various  strategies  for  producing 
computaole  formulations  of  D(X,I).  The  third  addresses  value 
standards  for  judging  a  computable  formulation  to  be  "good,'' 
or  for  judging  one  good  computable  formulation  "better"  than 
another.  And  the  fourth  issue  involves  the  selection  of  a 
sequential,  32-bit  word  length  environment  as  the  basis  for 
extrapolating  the  results  of  this  study  into  other  computing 
e.nv  i  ronments . 

In  view  of  anticipated  LCM  applications  in  m«ssomechanics , 
an  a  priori  judgment  of  "good"  will  be  assigned  to  a 
computable  formulation  of  DCN,I)  if  arithmetic  complexity 
analyses  predict  that  the  cpu  required  of  its  iaplmentat ion  in 
a  sequential  environment  is  independent  of  the  choice  of  s. 
We  see,  for  example,  chat  computation  of  the  sum  in  equation- 1 
via  a  ”00"  loop  statement  over  the  range  1  to  N  requires  cpu 
that  are  (linearly)  dependent  upon  the  site  of  M  and, 
therefore,  D^N.I)  is  not  a  "good"  computable  formulation  of 
D(M, I) .  Rationales  for  judging  a  given  good  formulation 
"better"  than  another  consist  of  appeals  to  issues 
articulating  memory  requirements,  parallelizeabilities ,  and 
the  range  of  (N,l)  pairs  over  which  the  given  formulation  is 
tractable  (i.e.,  avoids  overflow  and  is  otherwise  executable 
and  precise)  when  implemented  in  a  32-bit  word  length 
environment , 

Three  computable  formulations  are  produced  and  ranked 
according  to  these  issues  apparently  for  the  first  time  here. 

The  literatures  on  generalized  functions,  lattice  sums, 
et  al.,  offer  a  variety  of  strategies  for  otaining  computable 
formulations  of  0(N,I).  However,  there  seems  to  be  no 
previous  discussion  about  D(N,I)  with  respect  to  the 
computational  issues  expressed  above.  Consideration, 
therefore,  is  given  to  these  issues  as  they  relate  to  an 
immediately  perceivable  or  "naive"  strategy  for  producing 
computable  formulations  of  0(N,I).  Experiences  gained  from 
implementing  these  formulations  lead  to  development  of  an 
informed  strategy  for  producing  a  better  formulation,  and 
experiences  gained  thereupon  lead  in  turn  to  an  even  better 
f  ormulat  ion . 

Production  of  these  "better”  formulations  together  with 
the  rationale  given  for  judging  one  of  them  "best"  constitute 
the  main  points  of  this  paper.  The  significance  and  potential 
impact  of  the  paper  derive  from  the  proof  given  herein  that 
the  LCM  model  of  the  mesomechanics  problem  is  (1)  an  0(1)- 
system  of  equations  and,  therefore,  (2)  a  viable  alternative 
to  the  conventional  finite  element  model  (FEM)  because  the  FEM 
is  an  0(maxN)-system  where  0( I )  <<  O(maxN) . 
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8.  CONCLUSIONS  AND  RECOMMENDATIONS 

Conclusions 

We  have  presented  closed-form  approximations  for  the  frequency 
of  structural  systems,  both  with  and  without  damping.  As  such, 
these  solutions  are  valid  for  all  micromechanical  models  which 
determine  the  macromechanical  stiffnesses  Dy.  Contiguity  factors, 
fiber  shape  factors  and  misalignment  factors  only  affect  Dy,  and 
hence  are  implicitly  included  in  the  results  presented.  However, 
for  angle-ply  symmetric  laminates,  a  parametric  study  has  shown 
that  the  frequency  is  relatively  independent  of  the  micromechanical 
model  as  ? ong  as  classical  laminate  theory  applies  to  the  bonding 
between  layers.  It  is  further  pointed  out  that  FEM  models  do  not 
present  any  inherent  problem,  for  in  Chapter  3  a  solution  was 
obtained  for  the  general  discretized  problem. 

The  viscoelastic  analysis  presented  some  greater  than 
anticipated  difficulties.  Solutions  were  obtained  for  generic 
problems  in  which  the  damping  could  be  modelled  as  a  simple  damping 
matrix.  For  such  models,  however,  the  approximate  solutions 
presented  appear  very  promising. 

The  load  correction  method  has  been  extended  from  lattice 
elements  to  include  solid  elements  [1].  The  case  of  solid 
composite  elements,  however,  remains  to  be  demonstrated. 
Recommendations 

The  next  phase  for  the  LCM  approach  should  focus  on 
demonstrating  the  applicability  of  the  composite  solid  element  and 
the  completion  of  the  dynamic  analysis  using  this  element  to  obtain 
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natural  frequencies  for  composites.  This  could  be  accomplished 
using  either  the  methods  developed  by  Broome  in  Chapter  7  or  Reiss 
in  Chapter  3. 

Subsequent  suggested  work  on  sensitivity  analysis  could  focus 
on  various  simple  viscoelastic  models  of  orthtropic  materials  - 
including  Maxwell  and  Kelvin  models.  Once  done,  the  method  of 
Chapter  4  could  be  extended  to  these  viscoelastic  models. 
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